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1.0  Introduction 


A  physical  system  in  general  consists  of  hundreds  of  parts  which  are  connected 
in  various  forms.  The  governing  equations  for  the  system  can  be  formulated  in 
a  number  of  distinct  ways.  For  example,  Lagrange’s  approach[l],  which  is 
based  on  energy,  and  the  topology  approach [2],  which  is  based  on  casual 
patterns,  are  well  known.  The  governing  equations  obtained  are  nonlinear  in 
general  and  the  size  is  usually  very  large. 

Let  us  look  into  the  topological  approach.  The  development  of  governing 
equations  depends  on: 

(a)  the  structure  of  the  system  and 

(b)  the  properties  of  the  parts. 

The  laws  established  for  dealing  with  structures  are  Kirchoff’s  law, 
D’Alambert’s  principle,  etc.,  while  the  laws  for  treating  devices  are  Ohm’s 
law,  and  Hooke’s  law,  etc.  We  have  the  former  laws  of  constraints,  and  the 
latter  laws  of  performance. 

Nonlinear  phenomena  such  as  saturation,  dead-time  and  backlash[3]  frequently 
appear  in  the  laws  of  performance  and  nonlinear  junction  structure  or 
geometric  nonlinearity  quite  often  occur  in  the  laws  of  constra i nt[4] . 

After  the  governing  equations  has  been  developed,  the  first  treatment  to  these 
equations  is  usually  to  substitute  numerical  values  and  to  linearalize  them. 
The  solutions  of  course,  are  easy  to  obtain,  but  the  results  are  far  from  the 
reality  and  sometimes  lead  to  wrong  conclusions,  this  is  particularly  true  for 
most  practical  large  scale  systems [5] . 
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In  order  to  keep  the  original  modelling,  the  engineer  should  keep  symbolic 


analysis  as  long  as  he/she  can.  In  other  words,  to  substitute  numerical 
values  and  to  perform  linearizations  as  late  as  possible. 

However,  systems  to  be  investigated  become  larger  and  larger  and  computer’s 
capability  for  treating  symbols  has  not  been  improved.  Therefore,  there  is  a 
gap  between  the  system  modelling  and  qualitative  analysis.  there  is 
considerable  and  growing  interest [6]  in  establishing  symbolic  analysis 
especially  by  computer-aided  procedures  in  order  to  close  this  gap.  Recent 
examples  are  computerization  of  Mason’s  signal  flow  graph  and  Coates  flow 
formula  [7],  development  of  tree  enumeration  methods [8],  and  establishment  of 
parameter  extraction  techniques[9] .  However,  these  are  all  for  linear  models 
and  their  procedures  are  difficult  to  extend  to  the  nonlinear  domain.  It 
seems  that  there  is  a  missing  link  between  symbolic  analyses  of  linear  models 
and  nonlinear  models.  What  we  need  is  to  develop  a  general  scheme  to  cope 

with  the  fundamental  problem - symbolic  analysis  of  large  scale  linear  and 

nonlinear  systems. 


1.1  Whitehead’s  Universal  Algebra 

Since  Mason’s  time[7],  the  various  methods  developed  for  symbolic  analysis  are 
special  techniques  for  particular  types  of  questions  or  deliberately  coping 
with  certain  computing  machines.  No  general  foundation  has  ever  been 
established.  This  is  why  the  missing  link  from  the  nonlinear  models  to  linear 
ones  has  never  been  found. 
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Before  writing  the  well-known  classical  work,  Principia  of  Mathematica  with 
Russell,  Whitehead  wrote  a  book,  ’A  Treatise  on  Universal  Algebra’ [10]  which 
presents  a  thorough  investigation  of  the  various  systems  of  symbolic  reasoning 
allied  with  ordinary  algebra.  Whitehead  treated  three  basic  schools: 

(1)  Boolean  algebra, 

(2)  Hami Iton’s  formulation  and 

(3)  Grassmann  algebra. 

His  insight  was  far  ahead  of  his  time.  In  these  recent  years,  Boolean  algebra 
in  digital  computers  and  Hamilton’s  formulation  in  control  and  optimization 
theory  have  played  important  roles.  To  the  authors  of  this  report,  it  seems 

that  the  third  topic  of  Whitehead - Grassmann  algebra  is  a  powerful  tool  for 

investigating  systems  but  has  not  been  further  developed  since  his  time. 

In  a  preliminary  study,  the  authors  found  that  based  on  the  new  algebra  we  can 
cope  with  the  symbolic  investigation  of  large  scale  linear  systems  and  with 
further  and  deeper  studies,  we  may  approach  some  nonlinear  cases. 


1.2  Grassmann  Algebra  for  Determinants 

Let  us  briefly  introduce  some  of  Grassmann ’s [10]  basic  rules.  Suppose  we  want 
to  evaluate  the  determinant  of  a  matrix  M,  where 


ran  a12  a13  .  ?  nl 

I  a21  a22  a23  .  d2n  I 

I-  -I 

I-  I 

I-  I 

t-anl  an2  an3  .  ann-^ 


(1) 
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First  of  all,  each  column  is  expressed,  in  vector  notation,  a;,  or 

M  =  [  ai,  a2,  a3,  ,  an  ]  (la) 

where  a;  is  the  column  vector  of  the  i^  column.  The  determinant  of  M,  or 
A=det  M,  is  equal  to  the  exterior  product  of  the  vectors: 

A  =  ajAa2A...Aan  (2) 

where  A  is  a  symbol  denoting  exterior  of  outer  product  subject  to  the 
following  rules: 


ai 

>■ 

0) 

II 

O 

(3a) 

(ai 

;  A  >j)  .  (a,  A  a;)  a  0 

(3b) 

(a 

j  A  aj j  A  ak  =  a;  A  (aj  A  akJ 

(3c) 

ai 

A  (aj  ♦  ak)  =  (a;  A  aj)  ♦  (a;  A  ak) 

(3d) 

where  a,,  aj ,  and  ak  are  vectors  defined  in  (la). 

The  last  three  rules,  (3b)  to  (3d),  are  the  commutative  property,  associative 
property  and  the  distributive  property,  respectively,  while  the  first  rule, 
(3a)  is  unique. 

Let  us  use  an  example  to  illustrate  how  to  apply  the  Grassmann  rules  to  find 
the  value  of  the  following  determinant. 

Tan  a12  0  1 

Mi  =  |0  a22  a23|  (4) 

la31  a32  a33J 
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We  express  i ^  column  by  vector  a j ,  or 
al  =  (an ,  a31) 

a2  =  (a12<  a22»  a32)  (5) 

a3  =  (a23<  a33) 

It  is  noted  that  the  second  subscript  of  each  term  in  (5)  is  omitted.  Then 
the  determinant  of  is  evaluated  by  the  outer  product  of  the  column  vectors. 

Aj  =  det  Mi  =  ai  A  32  A  a3 

=  (1,3)  A  (1,  2,  3)  A  (2,  3) 

=  (11,  12,  13,  31,  32,  33)  A  (2,  3) 

=  (12,  13,  31,  32)  A  (2,  3) 

=  (122,  132,  312,  322,  123,  133,  313,  323) 

=  (132,  312,  123). 

(6) 

It  is  to  be  noted  that  when  we  performed  the  exterior  product  operations,  we 
have  used  rule  (3a).  Therefore,  11,  33,  122,  322,  133,  313,  323  become 

zero  in  various  steps.  Then  we  decode  equation  (6)  and  write  it  as: 

A1  *  (all  a32  a23»  a31  a12  a23»  all  a22  a33)  (7) 

the  second  subscript  of  each  element  in  (7)  is  added  according  to  the 
positions,  while  the  first  subscript  in  (7)  is  copied  from  (6). 

For  determining  the  sign  of  each  item  in  (7),  we  simply  check  whether  it  is 
odd  or  even  permutation.  For  example,  in  order  to  change  132  into  the  natural 
order  123,  we  need  one  interchange  of  the  positions  of  digits.  Therefore,  the 

term  132  should  be  with  a  negative  sign.  Similarly,  312  should  have  a  plus 
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sign,  because,  for  changing  it  into  a  natural  order,  we  have  to  interchange 
the  positions  of  digits  twice  which  means  even  permutation  and  a  positive  sign 
should  be  associated  with  the  item.  Finally,  we  have 


A1  =  ~all  a32  a23  *  a31  a12  a23  ♦  all  a22  a33  (8) 


1.3  State  Equation  Evaluation 


Consider  a  typical  state  space  model 


x  =  Ax  ♦  bu 
and  y  =  c^x 

where 

r  o  i  oi 

I  I 

A  =  |  0  0  1| 

I  I 

1-6  -11  -6J 

r  si 
I  I 

c  *  1 7 1 

I  I 

LOJ 


roi 
I  I 

b  =  joj  and 

I  I 

L1J 


(9) 

(9a) 


We  are  particularly  interested  in  the  evaluation  of  the  determinant 


1—4 

i/I 

-  A  ) 

(10) 

f  s 

-1 

01 

=  1  o 

s 

-1| 

l  6 

11 

s+6J 
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This  is  only  a  particular  case  of  (4),  if  we  interpret  the  i^h  column  as  a ; . 
Therefore, 

A  =  ai  A  a2  A  33  (4) 

=  -all  a32  a23  ♦  a31  a12  a23  ♦  all  a22  a33 
=  "(s)  (11)  (-1)  ♦  (6)  (-1)  (-1)  ♦  (s)  (s)  (s+6)  (4a) 

=  11s  +6  ♦  s^  ♦6s2 

=  s^  +6s2  ♦  Us  ♦  6  (11) 

which  is  what  we  expected. 

In  this  exercise,  the  Grassmann  algebra  enables  us  to  manipulate  a  practical 
symbolic  determinant  easily.  However,  when  the  system  is  large  and  when  we 
also  want  to  preserve  the  original  symbols,  the  problem  is  not  a  simple  one. 

The  evaluation  of  a  symbolic  determinant  can  be  extended  to  the  transfer 
function  evaluation.  Consider  the  following  identities 

r  P  Q  1 


det  | 

1 

=  det  P  det  (  T-RP"1!}  ) 

(12) 

L  R 

T  J 

r  p 

Q  1 

det  | 

1 

=  det  T  det  (  P-QT'iR  ) 

(13) 

L  R 

T  J 

These  two  identities  have  been  derived[ll]. 

Now  we  consider  the  set  of  state  equations  and  the  output  equation  in  (9). 
Our  problem  is  to  obtain  the  corresponding  transfer  funct  ion  Y(s)/U(s).  What 
we  need  is  to  interpret  the  two  identities  (12)  and  (13)  properly  and  apply 
grassmann ’s  algebra. 
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Let  R=-c^,  P=(sl-A),  Q=b  and  T=1 .  Substituting  the  new  notation  into  (12)  and 


(13)  we  obtain, 


1  ♦  c^(sI-A)~l  = 


dot  {(sI-A)+bcT) 
det  (sI-A) 


This  equation  was  derived  in  a  different  way  by  Patel  and  was  established  as 
above  and  used  in  the  deadbeat  systems  designed  by  Chen,  Tsay  and  Yates[ll]. 


The  left  hand  side  of  (14)  is  the  required  transfer  function  plus  one,  while 
the  right  hand  of  (14)  is  simply  a  ratio  of  two  determinants.  Therefore,  to 
evaluate  a  transfer  function  from  a  state  space  model  requires  the  evaluation 
of  two  symbolic  determinants. 


For  this  problem 


r  s  ■  -i  oi 

I  I 

det  (sI-A)  =  A<j  =  det  j  0  s  -1  j 

I  I 

l  6  11  s*6  J 


r  s  -1  0  1 

T  1  1 

det  { (sI-A) +bc '}  =  det  j  0  s  -1  j 

I  I 

L6+5  11+7  s+6  J 
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The  two  determinants  can  be  solved  by  the  Grassmann  algebra  method.  Therefore 


det  (sI-A)  =  -an  332  a23  ♦  831  a13  a23  ♦  an  a22  a33 
=  -(s)  (11)  (-1)  +  (6)  (-1)  (-1)  +  S2  (s+6) 

and 


det  {(sI-A)+bcT>  =  -aji  a32  a23  ♦  a3i  ai2  a23  ♦  a22  833 

=  -s(ll+7) (-1)  +  (6+5)  (-1) (-1)  +  s2(s+6) 

=  +  6s2  +  18s  +  11 

(17) 

Therefore,  the  desired  transfer  function  can  be  obtained  by  substituting  (16) 
and  (17)  into  (14). 

That  is, 


s^  ♦  6s2  ♦  18s  +  11 


1  +  cT(sI-A)’lb  =  - - - -  (18) 

s^  +  6s2  ♦  11s  +  6 

or 

.  7s  +  5 

J(sI-A)-lb  =  -- . . .  (19) 

s’  +  6s2  ♦  11s  +  6 


as  we  expected . 


1.4  Gain  Evaluation 


The  authors  of  this  report  generalized  the  utilization  of  Grassmann  algebra  to 
the  gain  evaluation  and  found  that  the  outgoing  branches  of  a  node  can  be 


interpreted  as  a  Grassmann  vector ,  and  the  determinant  of  a  cold  system  can  be 
obtained  by  taking  the  exterior  product  of  those  vectors. 

For  illustration,  consider  the  control  system  shown  in  Figure  1,  which  can  be 
represented  by  the  block  diagram  shown  in  Figure  2.  The  block  diagram  can  be 
expressed  by  a  matrix  equation: 


r-i  o  o 

0  0 

1  f  If  1 

1  I  1 

r-l/(  Rf  ♦  sLf  )  1 

1  1 

|Kg  -i  0 

0  0 

1  1  1 

1  1  Eg  1 

1  1  3  1 

1  1 

|  0  | 

1  f 

|  0  l/CRj+sLi)  -1 

0  -l/CR^sLx)!  |  Ia  |  = 

i  i  i 

1  0  1  Ef 

1  1 

|  0  0  Kt/(Js+b) 

-1  0 

1  1  Qm  1 

t  1  1 

1  1 

|  0  1 

1  1 

L  0  0  0 

Km  -1 

1  9  1 

J  L  ^  J 

1  1 
l  0  J 

(20) 


To  apply  the  Grassmann  algebra  to  equation  (20),  we  first  code  each  element  in 
the  matrix  equation  by  a  double  script  symbol  and  obtain 


r  an 

0 

0 

0 

0  1 
| 

r  if  l 

1  1 

r  h  1 

i  i 

1  a21 

a22 

0 

0 

1 

o  1 

1  i 
!  ^  ! 

i 

j  0 

1 

1 

1 

1  o 

a32 

a33 

0 

a35  1 
| 

1  Ia  1  = 

1  0 

1 

1  Ef 
| 

1  o 

0 

a43 

a44 

o  ! 

i 

1  Qrn  1 

1  | 

1  0 

1 

1 

1 

1 

L  0 

0 

0 

a54 

i 

a55  J 

1  1 
l  Em  J 

1 

l  0 

1 

J 

(20a) 

where  =  -1/  (  Rf  sLf  )  . 

We  can  now  draw  the  flow-graph  in  the  Coates  sense  which  is  shown  in  Figure  3, 
either  directly  from  the  block  diagram  or  from  the  matrix  equation. 

Let  us  consider  the  cold  system  first  which  means  that  the  system  is  without 
input.  Redraw  Figure  3  by  deleting  all  the  incoming  branches.  Then  we  have 
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Figure  4.  The  determinant  of  the  cold  system  can  be  found  from  the  exterior 
product  of  the  vectors  formed  by  the  outgoing  branches  of  each  node.  In  other 
words . 


Ac  s  (all/a2l)A(a22/a32)A(a33'a43))A(a44>a54)A(a35>a55)  (21) 


In  Grassmann’s  notation,  we  write 


Ac  =  (1»  2)  A  (2, 3)  A  (3, 4)  A  (4, 5)  A  (3, 5) 


Then  the  exterior  product  operations  are  performed  subject  to  the  Grassmann’s 
rules  shown  in  Section  1.2.  We  have 

Ac  *  (12 ,  13 , 23)  A  (3 , 4)  A  (4 , 5)  A  (3 , 5) 

=  (123, 124 ,134 ,234) A (4, 5) A (3, 5) 

=  (1234 , 1235 , 1245 , 1345 , 2345) A (3 , 5) 

=  (12345,12453)  (22) 


which  means 


all  a22  a33  a44  a55  ♦  all  a22  a43  a54  a35 

After  checking  even  or  odd  permutation,  the  sign  of  each  item  is  decoded. 
Therefore,  we  have 

Ac  =  ♦  all  a22  a33  a44  a55  ♦  all  a22  a43  a54  a35  (23) 

Substituting  the  original  values  gives 
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Fignro  1.  A  Control  Sjito* 


Block  Dli|na  of  tko  Control  Sjitu 

(I,-!.**.. 


-1 


kT 

(-1)  (-1)  (-1)  (-1)  (-1)  ♦  (-1)  (-1)  ( . )(Km)( 

Js«-b 

*T  Kffl 

-1 . 

(Js  ♦  b)  (R^  ♦  sLj) 


- ) 

Rl+sLi 


Which  is  the  denominator  of  the  gain. 


(24) 


If  we  are  interested  in  the  gain  between  say,  fim  to  Ef,  we  write  the  output 
equation  as  fol lows 


y  =  [0,  0,  0,  1,  0]  [If,  Eg,  Ia,  flm,  Ejh]^  (25) 

=  Jx 


By  us i ng  (14) ,  we  f i na I  I y  get  the  numerator  of  the  ga i n 

1  KT  -1 

N  =  -(Kg  x - x - x - x  (-1)  )  (26) 

Rl+sL^  Js  ♦  b  Rf+sLf 

Combining  (24)  and  (26)  yields  the  gain,  nm(s)/Ef (s) . 

As  we  have  seen,  the  new  method  based  on  Grassmann’s  algebra  is  not  only 
suitable  for  computers  to  deal  with  large  scale  systems,  but  also  preserves 
symbols  all  the  way.  While  Mason’s  formula  frequently  makes  mistakes  for 
large  scale  operations,  and  Coates’  formula  is  only  a  different  form  of 
Cramer’s  rule,  our  approach  is  more  accurate  than  that  of  Mason  and  more 
general  than  that  of  Coates. 

For  the  simplicity  of  generalization,  we  need  to  formalize  the  steps  from  (20) 
to  (26).  Since  any  single  input  and  single  output  (SISO)  signal  flow  graph  is 
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equivalent  to  a  linear  system  of  the  form 


(27) 

(28) 


Ax  =  bu 
y  =  cTx . 

where  A  is  the  Coates  graph  matrix,  u  the  input,  y  the  output,  b  the  column 
vector  consists  of  all  the  negative  gains  from  the  input  node  to  all  other 
nodes  in  the  flow  graph,  and  c^  is  the  row  vector  of  the  form 
[cj, . . . ,c j , . . . ,cn] ,  where  c;  is  the  gain  from  the  i^  state  variable  node  to 
the  output  node  y. 

We  have 

x  =  A_1bu, 
and  y  =  c^A^bu, 

the  transfer  function  or  the  gain  of  the  flow  graph 

G  =  cTA-lb  =  An/Ad.  (29) 

Compare  the  c^A"^b  with  the  c^(sI-A)"^b  in  (14),  we  have 

1  ♦  c^-lb  =  (  An  -  Ad  )/(  Ad  ) 
det  (A  +  bc^) 
det  A 

det  (A  bc^) 

or  1  ♦  G  =  .  (30) 

det  A 

For  elaboration,  examine  the  following  example. 
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Example  1 . 


Evaluate  the  gain  of  the  discrete  signal  flow  graph  shown  in  Figure  5. 

Using  matrix  notations  and  the  dummy  variables  aj  and  <12,  the  flow  graph  is 
equivalent  to  the  following  linear  system 


r  -1 

i 

-4 

0 1 

1 

r  Ai  i 

1  t 

f  -1  1 

1  1 

I  0  |  r (k) 

1  1 

i  z-1 

-1 

1 

0 1 

1 

11 

CM 

1  1 

ck 

1 

-1  j 

L  u(k)  J 

1  1 

L  0  J 

r  ai  1 
I  I 

u(k)  =[001]  |  a2  | 

I  I 

l  u(k)  J 


Then,  the  cold  system (flow  graph  without  the  input)  is  equivalent  to  the 
Coates  graph  shown  in  Figure  6. 

The  determinant  Ad  is  found  by  evaluating 

r  1  1  0  1 

I  I 

|  2  2  0  | 

I  I 

L  3  3  3  J 

symbol  ica  My 

Ad=(123)-(213) 

=  (-1) (-1) (-1)  -  (Z-l)(-dk)(-l) 

•  -(l.dkZ-l) 

To  found  the  numerator  An,  we  need  to  identify  the  b  and  cT  vectors 
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17 


18 


r  o  o  -l  1 


bcT  =  I  0  0  0  I 

I  I 

L  0  0  0  J 


Therefore, 


r  -1  -dk  -1  1 

T  I  ,  I 

A  -  bcT  =  |  Z-l  -1  0  | 

I  I 

II  ck  -1  J 

which  is  evaluated  symbolically  by  using  the  matrix 

r  1  1  1  1 

I  I 

|  2  2  0  | 

I  I 

L  3  3  3  J 
Therefore, 

An  ♦  Ad  =  (1  2  3)  -  (2  1  3)  ♦  (2  3  1)  -  (3  2  1) 

*  (-1) (-1) (-l)-(Z-l) (-dk) (-l)-(Z-l) (ek) (-l)-(l) (-1) (-1) 

•  -(UdkZ-l).(l.c|<Z-l), 

»n  *-(l*ckZ-l). 

From  (29), the  gain 

u(k)  An  l+c^Z”1 

G(Z-l)  = - =  —  = - - 

r(k)  Ad  l*dkZ”^ 

The  equations  (29)  and  (30)  will  be  generalized  in  the  next  chapter. 


19 


2.0  Symbolic  Analysis  of  Multivariable  Systems 


This  chapter  is  an  extension  to  the  last  one.  We  will  generalize  all  the  SISO 
results  to  Multi -inputs  and  Mu  I ti -outputs  (MIMO)  cases. 


2.1  Generalization  of  Crassmann  Algebra 


From  the  previous  Chapter,  for  a  SISO  state  equation,  we  have  the  result  (14) 


1  ♦  cT(sI-A)-1b 


det  {(sI-A) ♦be"'"} 


dot  (sI-A) 


For  an  n  inputs  and  m  outputs  system,  the  state  equations  are 


x  =  Ax  ♦  Bu 


Y  =  CTx 


The  transfer  function  matrix 


G(s)  =  CT(sI-A)-lB  =  [  gij(s)  ]  =  |  . 

I  • 


f  911 (s)  912 Cs)  •••  9lm(s)  1 
I  921 (s)  922 (s)  •••  92m (s)  I 


I 


where 


CT  = 


r  C1T  1 

I  c2T  | 


I 


I 


and 


t  9nl(s)  9n2(s)  •••  9nm(s)  J 


B  =  [bj,  b2,  -  •  • ,  bj  . 


(31) 


l  cnT  j 


Therefore, 

f  C1T(sI-A)-1  1 
|  c2T(sI-A)-l  | 

G(s)  =  !  ■  I  [  b!,  b2,  ...,  bm  ] 

1  '  T  ,  ! 

l  cnT(sI-A)-l  J 
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r  cJCsI-AHbi  ClT(sI-A)-lb2  • 
|  C2~(sl-A) C2^(sl-A) _^b2  . 


(sI-A)-Jh.  1 
(sI-A)-lbm  | 
I 


1  eJfsJ-AJ-lb!  : 


c„T(sI-A)-lb,,  J 


Then,  we  have 
1  ♦  9ij(s) 


1  ♦  Cj^"(sl-A)~^bj 

Aij/Ad» 


det  { (sI-A) +bjc  jT> 
det  (sI-A) 


where  Ajj  =  det  {(sI-A) *bj cj^}  -  Aj, 
and  Aj  =  det  (sI-A) . 


(32) 


(33a) 

(33b) 


Comparing  (32)  and  (33)  with  (14),  it  is  found  that  we  have  reduced  the 
evaluation  of  the  nXm  transfer  matrix  into  the  evaluation  of  nm  SISO  transfer 
functions. 


Example  2. 

To  illustrate  the  generalization  of  Grassmann  algebra,  consider  the  RLC 
network  shown  in  Figure  7.  Let  the  voltages  Vj  and  V2  be  the  inputs  and  the 
currents  Ij  and  I2  be  the  outputs.  The  loop  equations  are 

-Vi  ♦  IjR  ♦  (I1-I2)/sC  =  0, 
and  V2  +  I2(R+sL)  +  (l2-Ii)/sC  =  0. 
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rearranging,  we  have 


TR+ (1/sC) 

| 

-(1/sC)  1 

f  11  1  = 

ri  oi 

1  i 

r  Vi  l 

i  i 

L - (1/ sC) 

R+(l/sC)+sLj 

l  I2  J 

1  1 
10  -1J 

1  1 
L  V2  J 

r  yi  1 

n  oi 

f1*] 

Y  = 

1  1  = 

1  i 

' 

l  y2  J 

10  1J 

L  I2  J 

The  transfer  function  matrix  is  given  by 


«(•)  =  [  9ij  (S)  ] 


r  9n  (s) 
1  921 (s) 


912 (s)  1 


I 

922  (s)  J 


and  the  gjj(s)  is  given  by  (32) 


Ad  =  det 


rR+(l/sC)  -(1/sC)  1 


L-(l/sC)  R+(l/sC)+sLJ 


R+  (R+sL)  (1+RCs) 
sC 


_  r  i  i  r  i  o  l 

bici1  =|  I  [1  0]  =  |  | 

l  0  J  l  0  0  J 


Therefore, 


r  l+R+(l/sC)  -(1/sC)  1 

An  4  Ad  =  det  |  | 

L  -(1/sC)  R+(l/sC)+sL  J 


All 


(1+RCs) (s^LC)  -  1 

(sC)2 


911  =  An/Ad 

(1+RCs) (s^LC)  -  1 
(sC) [  R+(R+sL) (1+RCs)  ] 
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To  evaluate  the  g\2>  we  use 


r  o  i  r  o  oi 

b2C1T  =|  I  [1  0]  =  |  | 

t  -1  J  t  -1  0  J 


A12  +  Aj  =  det 


r  R+(l/sC)  -(1/sC)  1 

I  I 

L  -(!-*■  (1/sC) )  R+(l/sC)+sL  J 


A12  =  -(1/sC), 


Therefore 


912  =  A12/Ad  = 


-1 


R+(R+sL) (l+RCs) 


Similarly,  we  obta i n 

A21  =  1/sC» 

-(l+RCs) 


and 


A22  = 


The  transfer  matrix  is 


«(•)  = 


r  (l+RCs) (s^LC) -1  -1  1 

(sC) [R+(R+sL) (l+RCs)]  R+(R+sL) (l+RCs) 

1  -(l+RCs) 

L  R+(R+sL) (l+RCs)  R+(R+sL) (l+RCs)  J 


Examp  I e  3 . 

Consider  the  state  equations 
x  =  Ax  +Bu 
Y  =  CTx. 

where 
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r  i  i  0001 

|2  2  2  0  0j 
|  3  0  3  0  0  j 
|4  0  0  4  0| 
l  5  0  0  5  5  J 


and  the  determinant  is  given  symbolically  as 


(1  2  3  4  5) - (1  2  3  5  4)-(2  1  3  4  5)  +  (2  1  3  5  4)*(3  1  2  4  5)-(3  1  2  5  4) 


For  the  (sI-A+b2C2^) ,  we  must  use  the  symbolic  matrix 


r  1  1  0  1  0  1 

|2  2  2  0  0j 

|3  0  3  0  0| 

|4  0  0  4  4| 

l  5  0  0  5  5  J 
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and  the  determinant  is  given  symbolically  as 

(1  2  3  4  5)  —  (1  2  3  5  4)-(2  1  3  4  5) ♦  (2  1  3  5  4)* 

(3  1  2  4  5)-(3  1  2  5  4) - (4  2  3  1  5)  +  (5  2  3  1  4) 

For  illustration,  we  find  det(sI-A) . 

r  s+4  -i  o  o  oi 
|  5  s  -1  0  0  | 

(sI-A)  =  |  2  0  s  0  0  | 

|  1  0  0  s+4  -1  j 

L  2  0  0  4  s  J 

on  substitution, 

det(sI-A)  =  Ad  =  (s+4) (s) (s) (s+4) (s)-(s+4) (s) (s) (4) (-1)- 

(5)  (-1)  (s)  (s+4)  (s)  +  (5)  (-1)  (s)  (4)  (-1)* 

(2)  (-1)  (-1)  (s+4)  (s)-(2)  (-1)  (-1)  (4)  (-1) 
=  s5  ♦  8s4  ♦25s3  ♦  38s2  +28s  +8 
=  (s  ♦  l)2  (s  ♦  2) 3 

For  the  det{(sI-A)+b2C2^}, 

r  s+4  -1  0  1  0  1 

_  I  5  s  -1  0  0  | 

(sI-A)+b2C2'  =  j  2  0  s  0  0  j 

I  1  0  0  s+3  -1  j 

1  2  0  0  4  s  J 

on  substitution, 

det{ (sI-A) +b2C2T>  =  (s+4) (s) (s) (s+3) (s) - (s+4) (s) (s) (4) (-1) - 

(5)  (-1)  (s)  (s+3)  (s)  +  (5)  (-1)  (s)  (4)  (-1)  + 

(2)  (-1)  (-1)  (s+3)  (s)  -  (2)  (s)  (s)  (1)  (-1)- 
(1)  (s)  (s)  (1)  (s)  +  (2)  (s)  (s)  (1)  (-1) 

=  s3  ♦  7s4  ♦  20s3  ♦  31s2  ♦  26s  +8 
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&22  =  det{ (sI-A) ♦b2C2^}  -  det{(sI-A)} 

=  -s  (  s3  ♦  5s2  ♦  7s  ♦  2  ) 

=  -s  (  s  ♦  2  )  (  s2  ♦  3s  ♦  1  ) 

922 (s)  =  A22/Ad 

-(s) (s*2) (s2*3s*l)  -s(s2+3s+1) 

(s*1)2(s*2)3  (s*1)2(s*2)3 

using  the  different  (sI-A)+bjCj"^  matrices,  we  have 

det(sI-A*'  jc^)  =  s3  ♦  8s*  25s3  ♦  38s2  28s  ♦  8 

det(sI-A+b2CjT)  =  s3  ♦  9s*  ♦  29s3  ♦  42s2  ♦  28s  ♦  8 

det(sI-A+biC2^)  =  s3  ♦  7s*  ♦  21s3  ♦  33s2  ♦  26s  ♦  8 

Then 

An  =  0 

A12  “  s4  ♦  4s3  ♦  4s2  s  s2  (s  ♦  2) 

A21  *  -s(s3  ♦  4s2  ♦  5s  +  2)  a  -s(s*1)2(s*2) 
and  the  transfer  function  matrix  is  given  by 


r 


0 


G(s) 


(s*2)2 


s2  1 
(s+1)2(s+2)  I 
-s(s2*3s+1)  I 
(s*1)2(s+2)  J 


2.2  Multivariable  Signal  Flow  Graphs 


The  general  form  of  the  linear  equations  of  MIMO  signal  flow  graph  is 


Ax  =  BU 
Y  =  CTx. 
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where  A  is  the  Coates  graph  matrix,  U  the  input  matrix,  Y  the  outputs,  B  is  a 
matrix  with  the  i^  column  the  negative  gains  from  the  i^  input  to  all  other 
nodes  in  the  flow  graph,  and  is  a  matrix  with  the  row  the  gains  from 
the  state  variable  nodes  to  the  j^1  output  node  yj  of  Y. 

Which  implies  that  the  transfer  function  matrix 

G(s)  =  [  gjj(s)  ] 


F  ci^A~lbi  . . .  c]^A”lbm  1 

I  *1 

=  I  .1 

L  cn^A“^bi  . . .  cn"^A  *bm  J 

that  is 

gjj  =  Cj'A-ibj  =  Ajj/Arf  (34a) 

(  A-,:  ♦  A<j  )  det(  A  ♦  b;cjT  ) 

or,  1  ♦  gij  =  -  =  -  (34b) 

Ad  det  A 

Comparing  (34a)  and  (34b)  with  (29)  and  (30)  respectively,  we  will  find  that 
we  have  reduced  the  symbolic  gain  evaluation  of  an  m  inputs  and  n  outputs 
signal  flow  graph  into  the  gains  evaluation  of  nm  single  input  and  single 
output  signal  flow  graph. 

Example  4. 

Evaluate  the  transfer  function  matrix  of  the  discrete  signal  flow  graph  shown 
in  Figure  8. 


=  cta-ib 

r  CJ1  1 

I  •  I 


lcnTj 


A'1  [  blf  b2,  ...,  bm  ] 
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FIGURE 


The  2- inputs  and  2-outputs  signal  flow  graph  is  equivalent  to  the  following 
I i near  system 


-1  -Gi  -G2  0  0  0 

0 

01  f  X1  1 

r 

-1  0  1 

z-1  -1  0 

0  0  0 

0 

0  |  |  x2  | 

1 

0  0  | 

o  z-i  -1 

0  0  0 

0 

0  |  |  x3  | 

0  0  | 

0  0  Kx 

-110 

0 

0  |  |  x4  | 

1 

0  0  | 

0  0  0 

0  -1  -g3 

0 

0  1  j  x5  |  = 

1 

0  0  j 

0  0  0 

0  Z-l  -1 

0 

0  |  |  x6  | 

1 

0  -2  | 

0  0  0 

0  0  Z-l 

-1 

0  1  1  x7  | 

1 

-1  0  j 

0  0  0 

0  0  0 

k2 

-1  J  L  x8  J 

l 

0  0  J 

ui(k)  1 

r  o  o  o 

1 

1 

0  0  0  0  1 

u2(k)  J 

L  0  0  0 

0 

0  0  0  1  J 

r 

-1  -Gi  -G2 

0 

0  0  0  0 

1 

1 

l- 1  -1  0 

0 

0  0  0  0 

1 

1 

0  Z-l  -1 

0 

0  0  0  0 

1 

=  det  j 

0  0 

-1 

10  0  0 

1 

0  0  0 

0 

-1  -G3  0  0 

1 

1 

0  0  0 

0 

z-i  -10  0 

1 

1 

0  0  0 

0 

0  Z-l  -1  0 

1 

L 

0  0  0 

0 

0  0  K2  -1 

J 

cal ly, 

r 

1110 

0 

0  0  0  1 

1 

2  2  0  0 

0 

0  0  0  1 

l 

0  3  3  0 

0 

0  0  0  j 

1  =  det  | 

0  0  4  4 

4 

0  0  0| 

1 

0  0  0  0 

5 

5  0  0  1 

1 

0  0  0  0 

6 

6  0  0  j 

1 

0  0  0  0 

0 

7  7  0| 

L 

0  0  0  0 

0 

0  8  8  J 

=  (123 

45678)- 

(1 

2  3  4  6  5  7  8) 

-(2134 

(2  1  3 

4  6  5  7  8)'* 

(2  3  1  4  5  6  7  8) 

-(2314 

r  rx(k)  1 
I  I 

I  I 

t  r2(k)  J 


=  (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-l)-(-l) (-1) (-1)  (-1) (Z-l)  (-G3)  (-1) (-1)- 
(Z-1)  (-Gj)  (-1)  (-1)  (-1)  (-1)  (-1)  (-1MZ-1)  (-Gi)  (-1)  (-1)  (Z-l)  (-G3)  (-1)  (-1)* 
(Z-l) (Z-l) (-G2) (-1) (-1) (-1) (-1) (-l)^(Z-l) (Z-l) (-G2) (-1) (Z-l) (-G3) (-1) (-1) 
=  1  ♦  G3Z”1  ♦  GjZ-1  ♦  GiG3Z“2  ♦  G2Z“2  ♦  G2G3Z“3 
=  (1*GjZ“1+G2Z-2) (1*G3Z-1) 
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An  ♦  Ad  = 


r  -i  -Gi  -g2  -i  o  o  o  oi 

|  Z-1  -1  0  0  0  0  0  0  I 

|  o  Z”1  -1  0  0  0  0  0  j 

|  0  0  Ki  -1  1  0  0  0  I 

|0  0  0  0  -1  -G3  0  0  I 

|0  0  0  0  Z-1  -1  0  0  I 

jo  0  0  -1  0  Z_1  -1  0  I 

L  0  0  0  0  0  0  «2  -1  J 

Symbo I i ca I  I y , 


A11  +  Ad 


riii  iooooi 
|22000000| 
(033000001 
|00444000| 
(00005500  j 
|00006600| 
|00070770| 
L000000  8  8J 


=(1234567  8) — (1  234657  8) - (2  1345678)* 
(2  1  3  4  6  5  7  8)  +  (2  3  1  4  5  6  7  8)-(2  3  1  4  6  5  7  8)- 
(2  3  4  1  5  6  7  8)  +  (2  3  4  1  6  5  7  8) 

=  Ad  ♦  KiZ-2(1+G3Z-1) 

Therefore  An  =  KiZ-^l+GsZ"1) 

A11  K1Z'2 

911  = - = - : - " 

Ad  l+GiZ-i^Z"2 

Similarly,  we  obtain 

a12  =  K2(1*GiZ-UG2Z-2)  (UG3Z-1) 

A2i  =  2 (1+GiZ“1+G2Z"2) 
a22  =  (2K2Z-2)(UGiZ-Ug2Z-2) 

The  transfer  matrix 


r  KiZ*2  i 

I  . —  K2  I 

|  1*GiZ'1*G2Z-2  | 

!  o  1 

|  2  2K2Z-2  | 


L  1+G3Z~1  1+G3Z“^  J 
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Furthermore,  we  show  the  superposition  technique  for  tackling  MIMO  system 
problems  here.  Using  the  superposition  property  of  linear  systems,  different 
components  of  the  transfer  function  matrix  can  be  found  by  the  reduced  flow 
graphs  in  Figure  8. 

From  Figure  0(a),  the  gn(Z“l)  is  found  by  use 

f  —1  -Gj  -G2  0  1 

Ad  =  det  |  Z"1  -1  0  0  | 

|  0  Z"1  -1  0  | 

l  0  0  Kj  -1J 

and  is  evaluated  symbolically  by  using 

r  1  1  1  0  1 

122001 
I  0  3  3  0  1 

L  0  0  4  4  J 

Ij  s  (1  2  3  4) -(2  1  3  4)  ♦(2  3  1  4) 

=  (-1)  (-1)  (-1) (-l)-(Z-1) (-Cl) (-1) (-1) . (Z-l) (Z-l) (-C2)  (-1) 

=  1  ♦  CiZ"1  .  C2Z-2 

wi  th 

bT  =[-100  0]  and  cT  =  [0  0  0  1] 

An  *  Ad  's  evaluated  symbolically  by  using 

r  1  1  1  1  1 

|2  2  0  0| 

|  0  3  3  0  | 

l  0  0  4  4  J 

Au  *Ad  =(1  2  3  4)-(2  1  3  4) ♦  (2  3  1  4)-(2  3  4  1) 

=  (-1)  (-1)  (-1)  (-l)-(Z-l)  (— Gj)  (-1)  (-1MZ-1)  (Z-l)  (-G2)  (-1)- 
(Z-l)  (Z-l)  (KX)  (-1) 

=  Ad  ♦  KjZ-2 

A11  =  KlZ'2 

A11  K1Z-2 

gll(z-i)  -  — -  — — -  , 

Ad  UGjZ-l^Z-2 
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Similarly,  the  gj2  's  evaluated  by  using  a  2X2  matrix,  the  g2i  by  a  2X2  and 
the  g22  by  a  4x4.  Instead  of  using  an  8X8  matrix,  the  use  of  the 
superposition  property  of  linear  system  simplifies  the  gain  evaluation  to  a 
certain  extend.  This  method  has  some  similarity  with  the  aggregation  methods 
in  control  theories  and  worth  further  investigation. 

2.3  Flow  Graph  Decomposition 

When  the  number  of  nodes  in  a  flow  graph  is  large,  the  symbolic  evaluation  of 
the  determinant  of  the  graph  becomes  intractable.  A  recursive  decomposing 
technique  is  given  here. 

Given  any  matrix,  we  can  partition  it  into  the  following  matrix  block  form 

r  a  b  i 

I  I 

L  C  D  J 

where  A  is  partitioned  in  such  a  way  that  it  can  be  handled  by  the  Grassmann 
algebra  with  no  problems. 

Then,  as  (12), 

TAB! 

det  J  |  =  det  A  det(D-CA“^B) 

L  C  D  J 

Compare  the  matrix  CA“^B  with  the  C^(sI-A)“^B  in  (31),  we  have 

>  bm  ] 


cil  1 
c2T  1 


G  =  [  Sij  ]  =  CA-lB  =|  .  |  A-l  [  b1; 

J  if 


l  c  '  J 
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r  ciA*^bi  .  ciA-lbm  1 

_  |  cjA-lb}  .  .  .  C2A-lb:  I 

’  I  '  ,  •  I 

l  cnA  b^  .  cnA  bm  j 

That  is, 

det(  A+bicj^  ) 

1  ♦  9  i  j  = . £ -  (35) 

det  A 

(35)  implies  that  we  need  to  solve  a  MIMO  system.  In  terms  of  matrix  flow 
graph,  let 

r  a  b  i  r  Aii  a^2  i 

I  1  =  1  I 

L  C  D  J  L  A21  A22  j 

the  matrix  flow  graph  is  shown  in  Figure  10(a).  When  it  is  partitioned, 
Figure  10(b)  is  obtained.  As  far  as  the  determinants  concern,  we  can 
decompose  Figure  10(b)  into  Figure  10(c),  which  is  also  equivalent  to  Figure 
10(d). 

Since  we  partitioned  in  such  a  way  that  A^}  is  manageable,  there  is  no  problem 
in  evaluation  If  the  dimension  of  the  matrix  (A22“^21^11~^^12) 

is  still  too  large  to  manipulate,  we  can  repeat  the  algorithm  again. 

Examp  I e  3 . 

Consider  Figure  11(a),  evaluate  the  symbolic  gain  of  this  flow  graph. 

The  matrix  flow  graph  notation  is  shown  in  Figure  11(b)  and  the  Coates  graph 
matrix  is 
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FIGURE  11 
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To  find  the  determinant  Aj,  we  partition  the  matrix  into  the  form 

r  Ai  Bi  1 

I  I 

L  Ci  DX  J 

where 


r  -1  -G3 

0  2 

2ckl 

r  -1  -Gi 

-G2 

0  J 

1  z-i  -1 

0  0 

0 

i  z-i  -1 

0 

0  |  and  Di  * 

1  0  z-1 

-1  0 

0 

1  0  z-i 

-1 

0  1 

1  0  0 

k2  “1 

"dk 

1  0  0 

Kl 

-1  J 

L  0  0 

0  Z"1 

-1 

When  the  nodes  (xj,  X2,  x3,  X4)  are  grouped  together  as  one  matrix  node  and 
(x5<  x6»  x7»  x8*  x9)  as  another,  we  have  the  matrix  flow  graph  shown  in  Figure 
12(a).  As  far  as  the  determinants  concern,  we  can  decompose  the  matrix  flow 
graph  into  Figure  12(b).  Then 


r  Aj  b1  1 

Ad  =  det  |  | 

L  Ci  Dj  J 

=  det  Aj  det(  Dj  -  ) 

Since  =  [  0  ] 

Aj  =  det  Aj  det  Dj 
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The  determinants  are  evaluated  symbolically  by  using 


r  1  1  0  1  1  1 
r  1  1  1  o  1  |22000| 

|  2  2  0  0  |  and  |  0  3  3  0  0  | 
|0330|  |00444| 

L0044J  L00055J 


Ad  =  (  1+G1Z“^+G2z"2  )  [  2K2Z-2(l+ckZ-l)  -  (l+G3Z’l)  (UdkZ’l)  ] 


using  the  b  and  c^,  the  determinant  (  An  ♦  Ad  )  is  found  by  the  matrix 


I"  -1  -Gx  -G2  -1  0  0  0  0  0  1 

I  Z_1  -1  0000000  I 

|  0  Z-1  -1  0  0  0  0  0  0  I 

|  0  0  Kj  -1  1  0  0  0  0  I 

|  0  0  0  0  -1  -G3  0  2  2ck  I 

I  0  0  0  0  z-1  -1  0  0  0 

I  o  o  o-i  o  rl  -i  o  o  | 

|  0  0  0  0  0  0  K2  -1  -dk  j 

l  0  0  0  0  0  0  0  Z*1  -1  J 


Partition  into  the  form 


f  Ai*  Bi  1 

I  ,  ,  I 

L  Ci  Da  J 


where 


r  -l  -Gj  ~g2  -l  j 

Aa  =  |  Z'1  -1  0  0  |  and  Dx ’ 

I  o  z-i  -1  0  | 

1  0  0  Kx  -1  J 


r  Aj*  Bi’  1 

An  ♦  Ad  =  det  |  | 

l  Cx  Dx  J 

=  det  Aj*  det(  Dj  ’  -  C^Ai'-^*  ) 


Now,  det  Aa*  =  1  ♦  GiZ”1  ♦  G2Z-2  ♦ 


I  -1  -G3  0  2  2cki 

|  Z*-  -10  0  0 

=  |  0  z-i  -1  0  0 

|  0  0  K2  -1  -dk 

100  0  Z"1  -1  J 
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ciV-V  = 


|  0  0  0  0  0| 

I  931  0  0  0  0  I 

|  0  0  0  0  0| 

L  0  0  0  0  0  l 


Where 


931  = 


l+GiZ-1^"2 


i*g1z-1+g2z-2+k1z-2 


Therefore, 


Da*  -  cj  V-V 


r  -1  -G3  0  2  2ckl 

|  Z-1  -1  0  0  0  I 

|  -g31  Z-l  -10  0  I 

|  0  0  K2  -1  -dk| 

L  0  0  0  Z"1  -1  J 


det(Di  -C^  Aj  "^-Bi  )  =  2K2Z"2(l*ckZ~l)-(l+dkZ~l)  (1+G3Z“*)- 

931 (2K2) (l^ckZ-l) 

An  =  2K2Z-1(l*ckZ-1)[  K1Z-3-(l*G1Z-1*G2Z-2)  ]-K1Z”2(l^dkZ-1)  (l^Z"1) 


Then,  we  have 

u(k)  An  2K2Z-1(l+ckZ-l)  [K1Z-3-(UG1Z-UG2Z-2)]-K1Z“2(l+dkZ“1)  (l^Z"1) 
r  (k)  Ad  (l+G1Z-1*G2Z-2)  [2K2Z-2 (WckZ*i)  -  (UG3Z~i)  (l+dkZ'i)  ] 


To  illustrate  the  recursive  nature  of  the  algorithm,  we  apply  it  to  the 
evaluation  of  the  determinant  of  the  matrices  Dj  and  (Di-CjAi^Bi) .  Partition 
the  Dj  matrix  into 

r  a2  b2  i 

I  I 

i  c2  d2  j 


where 


r  -l  -g3  o  i 

1  i  1 

-  I  Z-l  -1  0  I 

1  1  1 
L  0  7-1  -l  I 


and  D2  = 


r  -1  -dkl 

1  1  1 

L  Z-l  -1  J 


det  =  det  A2  det(D2“C2A2”^B2) 
det  A2  =  -1  -  G3Z-1 

.  r  911  9121 

C2A2'l02  =  I  I 

L  0  0  J 


where 


911  = 


-(UG3Z-1) 


912  * 


2K2ckZ-2 

-(i*G3z-i) 


F  -l-gi 1  -dk "912  1 
det(D2-C2A2'*lB2)  *  det  |  | 

L  Z-l  -1  J 

s  (I*9ll)*.(dk+912)Z“1 

2K2Z-2(l+ckZ-l)  -  (1+G3Z_1) (l+dkZ_1) 
-(1*G3Z-1) 

det(A2)det(D2-C2A2~lB2)  =  2K2Z”2(i+Cj<2!“l)-(l+G3Z-^)  (l+dkZ~l) 

Partition  the  matrix  (Dj ’-Cj ’Ai ’)  into 


where 


r  A2’  Bo’  1 

i  ,  ,  I 

l  c2  d2  j 


r  -1  -g3  0  i 


A2’  =  |  Z*1  -1  0  |  and  D2*  = 

I  ,  I 

1-931  Z”1  J 


r  -1  -dki 

1  1  1 

L  Z"1  -1J 
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det  A2’  =  -1  -  G3Z-1 


C2,A2'-1B2> 


I”  911*  912*  1 

I  I 

10  0  J 


911 

912 


2K2(-g3i+Z-2) 

-(I+G3Z-1) 

2K2Ck(-931+Z-2) 

-(UG3Z-1) 


det(D2,-C2’A2,-lB2’) 


f  -1- 


det  | 


911 


L  Z-l 


-dk-912’  1 

I 

-1  J 


=  (1+911  )  +  (clk+9l2)^-1 

2K2Z_2(l+c|<Z‘1)-(l+d|<Z~l)  (l«fG3Z-1)-g3i  (2K2)  (l+ckZ-1) 
-(I+G3Z-1)  "" 


so. 


detA2,det(D2,-C2,A2’-lB2’)  =  2K2Z-2(l+ckZ-1)-(l+dkZ-1)  (1+G3Z~1)- 

931 (2K2) (l+ckz_1) 


43 


3.0  Symbolic  Gain  Evaluation  of  Hybrid  Systems 


In  the  previous  chapters,  <ve  have  addressed  the  symbolic  gain  evaluation  of 
analog  and  digital  systems  separately.  But  in  real  life,  it  is  the 
combination  of  both  analog  and  digital  systems  most  abundant.  The 
manipulation  of  this  hybrid  systems  are  more  complicate  and  special 
methodologies  are  needed.  For  example,  Sedlar  and  Bekey[ll]  proposed  a  direct 
method  which  allows  the  evaluation  of  input-output  relations  directly  from  the 
hybrid  signal  flow  graph.  Kuo  [12]  also  presented  a  modified  Mason’s 
formulation.  The  method  we  proposed  here  use  the  Grassmann  algebra  for  the 
symbolic  gain  evaluation  and  suitable  for  automated  implementations  in  digital 
computers . 

3.1  Fundamental  Difficulties 


Consider  the  hybrid  system  shown  in  Figure  13,  there  are  two  samplers  in  the 
position  marked  with  y^  and  y2  respectively.  To  evaluate  the  symbolic  gain  of 
this  kind  of  hybrid  systems,  it  is  found  that  the  evaluation  will  be  more 
simple  by  using  yj  and  y2  as  the  state  variables  of  the  flow  graph. 
Reconfigurate  the  flow  graph  shown  in  Figure  13  to  get  the  one  in  Figure  14, 


Then, 

y2  =  Gi  [  R+G3  yi*  ] 

(36) 

and 

yi  =  G2  Y2 

(37) 

where  yi=yi(s),  y2=y2(s)-  And  yi*j  y2*  are  t^ie  sampled  version  of  the  state 
variables  yj  and  y2  respectively. 


Therefore, 


yi/G2  =  Ci  [  P  ♦  C3  yi*  ] 
which  is  equivalent  to 

yi  =  g1g2  C  R  ♦  g3  yi*  ]  (38) 

Taking  samples  at  both  side,  (38)  becomes 

yi*  =  (GiG2R)*  ♦  (GiG2G3) *yi* 

where  (GiG2R)  *  =  [Gi  (s)G2(s)R  (s)]  *  and  (G^R)*  =  [Gi  (s)  G2  (s)  R  (S)  ]  * 
Note  also  that  (GjG2R)*  *  Gi*(s)G2*(s)R*(s) ,  in  general. 

Therefore 


(GiG2R)* 

**  =  (CjR)*  -  (GjG^* .  (40) 

1  -  (G!G2G3)* 
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Some  I i terature[12]  states  that  (40)  is  the  pulse  transfer  function  of  the 
hybrid  system.  But,  in  general, 

(GjR)*  t  Gi*R* 

and  (GjG^)*  t  (GiG2)*R*. 

which  implies  that  (39)  is  inseparable  into  the  form 

W*  =  (T*)R*  (41) 

The  inseparability  of  (40)  implies  that  the  flow  graph  shown  in  Figure  14 
cannot  be  replaced  by  an  equivalent  Z-transfer  function.  It  makes  the 
usefulness  of  (40)  very  limited. 


To  alleviate  this  situation,  we  can  add  an  input  sampler  to  the  flow  graph,  it 
is  shown  in  Figure  15(a).  Then,  the  state  equations  becomes 


,2  »  Cl  [  R*  .  G3  ,i*  ] 
and  yi  =  G2  y2 

proceed  as  mentioned,  we  have 


W*  =  Gi*R*  ♦  (GjG3) * 

=  [Gi*  -  (GiG3)*  ■ 
=  T* (s)  R* 


(GxG2)*R* 

1  -  (G!G2G3)* 

(G1C2)* 


1  -  (GjG2G3) * 


]  R* 


using  (44),  the  simplified  flow  graph  with 


(G1G2)*(GiG3)* 

T*(»)  =  GX*  ♦  - 

1  -  (G!G2G3)* 


(42) 

(43) 


(44) 
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is  obtained  and  shown  in  Figure  15(b).  To  get  the  Z-transfer  function  of  the 
system,  use  the  following  substitution 


T(Z)  =  T*(s)  I s=  ( | nZ) /T 

Hereafter,  for  the  reason  of  implementations,  we  assume  that  there  are  always 
samplers  at  the  inputs  of  the  hybrid  systems  such  that  the  Z-transfer  function 
of  the  system  can  be  easily  obtained. 


3.2  Symbolic  Evaluation  Algorithm 


Consider  the  example  in  section  3.1,  the  state  equations  (42)  and  (43) 


y2  =  Ci  C  R*  -  C3  yi*  ] 

and  y^  =  G2  y2 

can  be  written  as  the  following  matrices  form 

r  yi  1  r  0  g2  i  r  yi  1  r  o  o  1  r  yi*  1  r  o  i 

11  =  1  II  I  ♦  I  II  I  ♦  I  l 

l  y2  J  L  0  0  J  L  y2  J  l  G^  0  J  l  y2*  J  L  Gj  J 

F  yj*  1 

w*  =  [  0  1  3  I  I 
t  y2*  J 


For  a  hybrid  system  with  n  samplers,  if  we  consider  the  input  of  the  i^*1 
sampler  as  y;(s),  the  i^  state  variable,  and  the  output  of  the  sampler  as 
y i * (s) ,  the  generalised  form  of  the  above  equation  is 
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(45) 

(46) 


y(s)  =  A1(s)y(s)  ♦  A2(s)y*(s)  ♦  b(s)R*(s) 

w*(s)  =  c*Ty(s) 

where  y(s)  =  [y^s),  y2(s),  yn(s)]T; 

y*(s)  =  [yi*(s),  y2*(s),  yn*(s)]T; 

Aj  (s)  and  A2(s)  are  nXn  matrices  and  b(s)  is  nXl; 

R*(s)  and  w*(s)  are  the  sampled  input  and  output; 
c*T(s)  =  C2*(s) . 

(45)  is  equivalent  to 

[I-Ai(s)]y(s)  =  A2(s)  y*(s)  ♦  b(s)R(s) 
or 

y  (s)  =  [I-AX  (s)]-l  A2(s)  y*  (s)  *  [I-AjCs)]-!  b(s)  R(s) 
=  Gn(s)  y*  (s)  ♦  G12(s)  R*(s) 

where  Gn(s)  =  [I-AjCs)]-1  A2(s); 

Gi2(s)  =  [I-AjCs)]-!  b(s). 

Sampling  the  both  sides  of  (47),  we  have 

y‘(s)  =  Gu*(s)  y* (s)  ♦  G12*(s)  R* (s) 

which  is  equivalent  to 

[I-Gn*(s)]y*(s)  =  G12* (s)  R* (s) 
or  y*  (s)  =  [1-Gn‘Cs)]-1  G12*  (s)  R*(s) 


(47) 


51 


Therefore,  from  (46) 


w*(s)  =  c*T(s)  [I-G11*(s)]“1  G12*(s)  R*(s) 

S  i  nee 

F(Z)  =  F*(s) |s=(|n2)/T, 
substituting  s=(lnZ)/T  into  (48),  we  have 

,*(Z)  =  c*T(Z)  <Si2*CZ)  R'CZ) 


r  GiG2  i 

G12 (s)  =  [  I-Al(s)  ]"lb(s)  =  |  I 

L  Gj  J 


Therefore, 


r  (GjG2) *1 
Gl2*(s)  =  |  | 

l  (Gl)*  J 


To  evaluate  the  gain  symbolically,  use  (51)  with 


cT(Z)  -  [  o,  1  3, 


r  GiGsCZ)! 

b(Z)  =  G12(Z)  =  |  | 

l  Gi(Z)  J 


=  det  [  I-GU(Z)  ] 

T  l-GiGsCsCZ) 

=  det  | 

l  -G!G2(Z) 

=  1  -  G1G2G3(Z) 


An  =  det  X  I-Gn(Z)*bcT  ] 

r  i-g1g2g3 (Z)  G!G2(Z)i 

=  det  |  I  -  Ad 

l  -GjG^Z)  UQi(Z)  J 


=  [  l-G^CZ)  ]Gi(Z)  ♦  G1G2(Z)G1G3(Z) 


Therefore, 


G(Z)  =  An/Ad 


=  Gl (Z)  ♦ 


GiG2(Z)G1G3(Z) 
1  -  GjG2G3(Z) 
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We  can  summarized  the  algorithm  as: 

1.  For  each  sampler  in  the  signal  flow  graph,  assign  a  state  variable;  the 
input  of  the  sampler  is  a  continuous  variable  y; (s),  while  the  output  of  the 
sampler  is  the  discrete  variable  y;*(s).  Notice  that  we  assume  that  there  is 
always  an  input  sampler  and  we  need  not  assign  a  state  variable  to  the 
sampler. 

2.  Reconfigurate  the  signal  flow  graph  such  that  there  are  only  ingoing  and 
outgoing  branches  from  the  state  variables  defined  in  (1). 


3.  the  state  equation  of  the  resultant  flow  graph  is  of  the  form 


r  yi  1  r  Anl(s)  A12J(s) 

j  y2  j  j  A2i1(s)  A221(s) 


••  AlnJ(s)  i  r  yi  i 

• •  A2nX (s)  I  I  Y2  I 

•  II*  I 

I  I  •  I 


L  yn  -1  l  Anl*(s)  An2^(s)  •••  Anr^(s)  ^  ^  yn  J 


r  AiX2(s)  A122(s)  . 
I  A212(S)  A222(s)  . 


I  Anl2(s)  An22(s)  . 


••  Ain2(s)  1  r  y*i  1  r  b2 (s)  1 
••  A2n2(s)  |  |  y*2  I  I  b2(s)  | 


•  .  I  I  •  I  I  ■  I 

••  Ann2(s)  J  L  y*n  J  L  bn(s)  J 


I  r*(s) 


*(s)  =  [  cjCs),  c2(s),  . 


r  yi*  i 
I  y2*  I 

,  C„(s)  ]  I  •  1 

I  .  I 
I  •  I 
I  yn*  J 


If  the  i j ^  elements  of  the  matrices  A*  (s)  and  A2(s)  be  Ajjl(s)  and  Ajj2(s), 
then  Ajjl(s)  is  the  weight  of  the  branch  going  from  the  input  of  the  j^h 
sampler  to  the  input  of  the  i^  sampler.  A;j2(s)  is  the  weight  of  the  branch 
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going  from  the  output  of  the  sampler  to  the  input  of  the  i^b  sampler. 

And,  bj (s)  is  the  weight  of  the  branch  going  from  the  output  of  the  input 
sampler  to  the  i^b  state  variable  node.  Cj(s)  is  the  weight  of  the  branch 

going  from  the  the  output  of  the  i^b  sampler  to  the  output  node  w(s) 

4.  Let  Gn(s)  =  [  I-Al(s)  ]"*  a2(s), 

and  G12(s)  =  [  I-Al(s)  ]'l  b(s). 

then  the  output 

w*(s)  =  c*T(s)  [  I-Gn*(s)  ]-l  G12*(s)  r* (s) 

and  the  Z-transfer  function 

G(Z)  =  cT(Z)  [  I-Cn(Z)  ]-l  G12(Z). 

5.  To  evaluate  the  Z-transfer  function  symbolically,  we  can  use  the  Grassmann 
algebra  defined  in  chapter  1  and  the  identity  (51) 

det[  I-Gn(Z)+G12(Z)cT(Z)  ] 

1  ♦  G(Z)  = . — - . — 

det  [  I-GU(Z)  ] 
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4.0  Decoupling  of  Multivariable  Systems 

In  this  chapter,  we  will  explore  the  problem  of  decoupling  a  MIMO  system.  We 
will  also  concentrate  in  systems  described  by  square  transfer  matrices  because 
any  system  described  by  non-square  transfer  matrices  can  be  decomposed  into 
several  systems  described  by  square  matrices.  The  main  reason  to  the 
discussion  of  decoupling  stems  from  our  concept  toward  MIMO  systems  design. 
We  want  to  decouple  the  MIMO  system  first  so  as  to  obtain  an  one- to  one 
correspondence  among  the  inputs  and  the  outputs.  Once  the  relation  is 
established,  the  well  developed  SISO  design  methodologies  can  then  be  applied 
to  the  systems. 

4.1  P-Constra i ned  and  V-Constra i ned  System  Model 

In  the  previous  Chapters,  we  use  Grassmann  algebra  to  obtain  the  following 
transfer  matrix, 

Y(s)  =  G(s)  U(s)  (52) 


If 


r  911  (s)  912(S)  •••  9ln(s)  1 
I  921  (s)  922 (s)  •••  92n (s)  I 
G(s)  =  |  .  .  .| 

I  •  -I 

l  9nl (s)  9n2(s)  •••  9nn(s)  J 

then  the  system  is  equivalent  to  the  block  diagram  shown  in  Figure  16. 
According  to  Mesarov ic[13] ,  the  structure  shown  in  Figure  16  is  known  as  a  P- 
constrained  system.  For  the  simplicity  of  further  explanation,  we  change  the 
notations  in  (52)  to  that  in  (53) 
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FIGURE  16 


CP(s)  =  P(s)  Mp(s) 


(53) 


where  Cp(s)  is  the  MIMO  output  matrix;  Mp(s)  the  input  matrix  and  P(s),  with 
P(s)=G(s),  as  the  P-constrained  system  matrix. 


There  is  another  basic  structure  of  MIMO  systems  known  as  V-constra i ned 
systems.  The  block  diagram  of  this  kind  of  systems  is  given  in  Figure  17(a). 
The  main  reason  to  investigate  the  P-  and  V-  constrained  system  structures  is 
that  they  need  different  decoupling  strategies  to  obtain  effectiveness  and 
simplicities.  All  these  will  be  explored  in  the  later  sections.  If  we  let, 
as  Figure  17(b) , 


r  Vn  C  0 

|  0  V22  0 

H  =  |  .  0  V33 

L  0  0  6 


and  K  = 


nn 


r  0  V12  V13 

I  V21  0  v23 

|  V31  v32  0 

i  vni  vn2  : 


•  Vm  1 

•  v2n  | 

•  v3n  I 

•  I 

.  0  J 


then,  we  have 

Cy  =  H  My  +  H  K  Cy 
or  CV  =  [  I  -  H  K  ]-l  H  My 


(54) 

(55) 


Comparing  (53)  and  (55),  if  both  systems  are  subjected  to  the  same  input 
signals,  i .e,  Mp  =  My,  then  the  condition  to  obtain  identical  responses  will 

be 


P  =  [  I-  H  K  J*1  H 


(56) 
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If  the  matrices  P  and  H  in  (56)  are  nonsingular,  then  taking  the  inverse  of 
(56) ,  we  have 

P-1  =  H-l  [  I-HK  ] 

=  H-l  -  K  (57) 

Since 

,  adJ  P 

P-1  - - , 

det  P 

where  adj  P  is  the  adjacent  matrix  of  P. 

Let 

r  Au  A12  ...  Aln  1 

*21*22  •••  A2n  | 

adj  P  =  |  .  .  | 

I  •  I 

l  Anl  An2  • - ■  Ann  J 

then  (57)  implies 


1 

f  Aai  o  o  . . . 

|  0  A22  0 

i 

0  1 

o  I 

1 

det  P 

•  o 

•  o 

•  o 

1 

Ann 

(58a) 

-1 

r  0  A12  A13  • 

|  A2i  0  A23  . 

|  .  .  . 

..  Aln  1 

••  A2n  1 

1 

det  P 

L  Anl  An2  • 

•  •  Ann 

(58b) 

If  we  denote  the  diagonal  matrix  constructed  by  taking  the  diagonal  elements 
of  any  matrix,  P,  as  [P] q ,  then  (58a)  and  (58b)  are  simplified  as 


H-l  .  [  P-1  ]D 
and  K  =  (  P-1  ]D  -  P-1 
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f59a) 

(59b) 


From  (56),  if  all  the  matrices  are 
whether  the  MIMO  system  is  P-  or  V- 
matr i x  P,  there  are  always  matrices 
vice  verse.  But  we  are  going  to 
that  whether  a  system  is  P-  or  V- 
structure  of  the  system. 


invertible,  it  is  impossible  to  identify 
constrained.  It  is  because  for  any  given 
H  and  K  such  that  (56)  is  satisfied,  or 
show,  through  the  following  two  examples, 
constrained  is  determined  by  the  physical 


Example  4 

Examine  Figure  18,  it  is  the  schematic  of  a  fluid  mixing  control  process.  The 
same  fluid  with  different  temperature  Tj  and  T2  are  inputed  from  A  and  B 
respectively.  It  is  hoped  that  the  fluid  with  temperature  T  and  flow  Q  is 
outputed  from  C.  The  problem  is  to  obtained  a  preseted  T0  and  Q0  by 
controlling  the  flows  and  Q2  through  the  valves  in  1  and  2. 

The  flow  of  the  fluid  is  preserved,  that  is 

Q  =  Ql  +  Q2 

In  steady  state,  let 

Ho  =  Hio  +  H20 

where  Qj0  and  Q2o  are  ^e  steady  state  valves  of  and  Q2  respectively. 

The  steady  state  error 


Qlo  &Hl  Q20  &H2 

•1  =  AQ/Qo  =  - . —  ♦ . -  (60) 

Ho  Hio  Ho  H20 
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where  AQ,  AQi  and  AQ2  are  the  deviations  from  Q0,  and  Q20  respectively 


Qlo  AQl  Q20  AQ2 

Ki  - - ,  K2  = - i  qi  = - ,  q2  = - • 

Qo  Qio  Qo  Q20 

Since  Kj  +  K2  =  1, 


then 

ei  =  Ki  qi  +  (i-Ki)  q2  (61) 


where  q^  and  q2  are  controlled  by  regulators  Rj  and  R2  respectively. 


From  thermo  equilibrium  considerations,  we  have 

(Ql  *  Q2K  =  QlTl  +  Q2T2 

Since  T  =  T(Q1(  Q2) 


AT  =  (T  -  T0) 

3T  9T 

= - AQj  + - AQ2 

0Ql  3Q2 


QloQ2o 

'oj' 


aq2 

(t2  -  Tl)  [  - 

Q20 


AQ2 


fiQi 

. ] 

Qlo 


where  Tm  = 


QloQ2o 

...... 


(t2  -  Tl)  • 


AQi 

Qlo 


Then  the  steady  state  temperature  error  is 


AT 


®2  = 


Tm 


AQ2 

Q2o 


AQi 

Qlo 


=  q2  -  qi 


(63) 
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From  (61)  and  (63),  we  obtain  Figure  19  and  the  following  matrix  equation 

r  1  1  f  Ki  l-Ki  1  f  qj  1 

11=1  III 

L  62  J  11  -1  J  L  q2  J 

which  is  a  P-constra i ned  system  structure. 

Examp  I e  5 

Examine  Figure  20,  the  objective  is  to  control  the  high  of  the  fluid,  h,  and 
the  output  flow,  Qq,  by  varying  yj  and  y0;  where  y;  and  yG  are  quantities 
proportional  to  the  variable  cross  areas  of  the  input  and  output  valves 
respectively.  The  input  flow  of  the  fluid  is  denoted  by  Qj  and  the  surface 
area  of  the  container  be  a  constant  A. 

We  have 

AQ o  =  C11  Ay,  ♦  C12  Ah  (64) 

where  Cji  and  C]^  are  positive  proportional  constants. 

S  i  nee 


A  dh  =  Qj  dt  -  Qq  dt 


After  a  time  t  elapsed, 


1  ft 

h  =  —  |  (Qj  -  Q^dt  ♦ho 

A  J  o 


where  hQ  is  the  initial  high  of  the  fluid. 
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That  is, 

1  ft 

Ah  =  h  -  h0  =  —  I  (AQ;  -  AQ0)dt 
A  Jo 

by  AQ j  =  Cj  Ay; 

where  Cj  is  a  positive  proportional  constant. 

We  have 

1  ft 

Ah  =  —  I  (C  j  Ay  j  -  AQ0)dt  (65) 

A  Jo 

Taking  the  Laplace  transforms  of  (64)  and  (65),  we  obtain 

AQq(s)  =  Cn  [  Ay0(s)  ♦  (  Cl2/Cn  )Ah(s)  ] 

Ah (s)  =  (C ; /As)  [  Ayj(s)  -  (l/Ci)AQ0(s)  ] 

which  is  equivalent  to  the  following  matrix  equation 

r  aQq  i  r  Ci i  o  1  r  Ay0  i  r  c^  o  i  r  o  (Ci2/Cii)i  r  aQq  i 

I  I  =  I  I  I  I  ♦  I  II  II  I 

L  Ah  J  10  (Cj/As)  J  L  Ay;  J  L  0  (C , /As) J  L(-1/Cj)  0  J  L  AQq  J 


Comparing  with  (54),  we  have 


(Cl2/Cll)  1 


f  Cjj  0  1  TO  (C12 

H  =  |  | ,  and  K  =  | 

l  0  (Cj/As)  J  l  (-1/Cj)  0 


That  is,  it  has  a  V-constra i ned  system  structure.  The  block  diagram  of  this 


system  is  shown  in  Figure  21. 
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FIGURE  20 
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4.2  Symbolic  Conversions 


When  we  need  to  convert  a  P-constrained  mathematical  model  to  a  V-constra i ned 
one,  we  can  use  the  Grassmann  algebra  so  that  symbolic  evaluation  is 
accomplished.  By  examinating  (59a)  and  (59b),  once  we  obtain  the  P-*,  the  H-* 
and  K  are  easily  constructed.  Therefore,  the  symbolic  evaluation  becomes  a 
symbolic  matrix  inverse  procedure. 

Examp  I s  5 
Consider  C  =  P  M, 


where 


r  i 

(s-l) 


P  = 


1 

(s+4) 

1 

(s+2) 


1 

L  (s+4) 


(s+3)  J 


The  block  diagram  of  this  P-constrained  structure  is  shown  in  Figure  22(a). 
To  obtain  the  H  and  K  used  in  the  V-constra i ned  description  (54),  we  need  to 
symbolically  evaluate  the  inverse  of  P. 

Using  Grassmann  algebra  with  =  B  =  I,  that  is 

C1T  =  [10  0],  c2T  =[010],  c3T  =  [001], 
and  bj  =  [1  0  0]T,  b2  =[0  1  0]T,  b3  =  [0  0  1]T. 
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then 


1 

Aj  =  det  P  =  - 

(s+1)  (s+2)  (s+3) 

and 

r  -(s+i)  (s+2)  i 

I  (s+l)  -  0  I 

i  (s+4)  | 

P-1  =  |  0  (s+2)  0  | 

I  I 

|  -(s+1) (s+3)  (s+1) (s+2) (s+3)  | 

I  . . :  •  (s+3)  I 

L  (s+4)  (s+4)2  J 


then  from  (59a)  and  (59b) ,  we  have 

r  (s+i)  o  o  i 

i  1  1 

H-1  =  |  0  (s+2)  0  | 

I  I 

l  0  0  (s+3)  J 


r  (s+i) (s+2)  i 

I  o  -  0  I 

I  (s+4)  j 

I  I 

K  =  |  0  0  0  i 

|  (s+1) (s+3)  -(s+1) (s+2) (s+3)  ! 

|  . —  0  j 

L  (s+4)  (s+4) 2  J 


The  resultant  block  diagram  of  the  equivalent  V-constra i ned  system  is  given  in 
Figure  22(b) . 


To  convert  Trom  tne  V-constra i ned  structure  into  the  P-constra i ned  one,  we 
can,  of  course,  use  (56).  But  it  is  simpler  to  start  from  (57).  It  is 
illustrated  by  the  following  example. 
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Example  7 


Consider  the  following  V-constra i ned  structure  with 


r  1 


H  = 


(*♦2) 

0 


(»♦*) 


and 


K  = 


0 

0 


(s«>6)  J 

-1  1 
(s+5) 

0 


0 
-1 

l  (s>3)  J 

The  block  diagram  of  this  system  is  given  in  Figure  23(a).  By  (57) 


P-1  =  H-l  -  K  = 


(»*2) 


1  1 

0  - 

(s+3) 

(s+4)  0 

1 


(*♦3) 


(s+6) 


J 


Once  again,  we  can  use  the  Grassmann  algebra  to  obtain  the  P  symbolically. 
Proceeds  as  the  previous  example,  we  have 
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r  i  1  -11 

I -  - - - - . .  i 

j  (s+2)  (s+2)  (s+3)  (s+4)  (s+5)  (s+6)  (s+2)  (s+5)  (s+4)  j 

I  I 

j  1  j 

P  =  |  0  -  0  | 

I  (s+4)  | 

I  I 

|  0  -1  1  | 

I  ! 

L  (s+3)  (s+4)  (s+6)  (s+6)  J 


The  resultant  P-constra i ned  structure  is  shown  in  figure  23(b). 

4.3  P-Constra i ned  and  V-Constra i ned  Decoupling  Structures 

Consider  (53)  again 

Cp(s)  =  P(s)  Mp(s)  (53) 

When  the  elements  pjj(s)=0  for  all  i*j,  that  is,  P(s)  is  a  diagonal  matrix,  we 
say  that  the  system  (53)  is  decoupled.  For  an  open  loop  system,  decoupling 
can  be  easily  achieved.  For  example,  examine  Figure  24(a),  it  is  a  P- 
constrained  system.  If  we  add  a  feedforward  path  in  this  open  loop  structure, 
as  shown  in  Figure  24(b),  the  system  will  be  decoupled  and  is  shown  in  Figure 
24(c). 

For  an  open  loop  V-constra i ned  system,  as  the  one  given  in  Figure  25(a),  we 
can  add  a  feedback  path  to  cancel  out  the  K  matrix  as  given  in  (54).  The 
feedback  path  is  shown  in  Figure  25(b)  and  the  decoupled  system  is  shown  in 
Figure  25(c) . 

But, in  practice,  we  use  a  lots  of  closed  loop  systems.  The  presentation  of  a 
simple  closed- loop  system  decoupling  method  is  the  objective  of  this  section. 
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FIGURE  25 
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To  decouple  the  closed  loop  system  shown  in  Figure  26(a),  we  can  put  a 
decoupling  structure  N  between  the  compensator  R  and  the  plant  P,  as  shown  in 
Figure  26(b).  Notice  that  the  matrices  R  and  F  are  unknown  yet,  but  we  can 
impose  the  condition  that  R  and  F  are  diagonal.  It  is  because  once  the  loop 
is  decoupled,  we  need  only  to  determine  the  corresponding  diagonal  elements  in 
R  and  F  by  SISO  design  methods. 

Since  the  couplings  of  MIMO  systems  can  be  divided  into  P-  or  V-  constrained 
structures,  as  illustrated  in  the  last  sections,  we  can  use  the  same 
structures  in  the  decoupler  N.  That  is,  we  assume  that  there  are  P- 
constrained  ana  V-constra i ned  decoupler  structures.  These  structures  of 
decouplers  are  illustrated  in  Figure  27(a)  and  (b) . 

Now,  by  the  assumption  that  R  and  F  are  diagonal,  we  want  to  obtain  the 
condition  that  N  will  decoupled  the  closed  loop  system  shown  in  26(a).  From 
the  figure,  we  have 

C  =  PNR(M  -  FC) 

which  is  equivalent  to 

[  I+PNRF  ]C  =  PNRM 
or,  C  =  [  I+PNRF  J"1  PNRM 
=  G  M 


where  G  is  the  transfer  matrix  of  the  closed  loop  system. 
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That  is, 


G  =  [  I+PNRF  ]-!  PNR  (66) 

We  want  the  closed  loop  system  decoupled,  which  is  equivalent  to  require  that 
G  is  diagonal.  Since  R  and  F  are  diagonal,  we  want  to  prove  that  G  is 
diagonal  if  and  only  if  the  matrix  PN  is  diagonal. 

Firstly,  we  want  to  prove  that  if  PN  is  diagonal  then  G  is  diagonal. 

since  the  addition  and  multiplication  of  diagonal  matrices  are  also  diagonal, 
the  matrices  W  =  I+PNRF  and  Z  =  PNR  are  diagonal.  Then,  (66)  becomes 

G  =  W-l  Z 

Since 

[  wno  o  ...  o  l*1  r  i/wu  o  o  ...  o  i 

|  0  W22  0  .  .  .  0  I  I  0  l/w22  0  ...  0  I 

w-1  =  I  •  •  .  |  =  I  .  .1 

I  ■  -I  I  •  -I 

10  0  ...  wnnJ  L  0  0  ...  l/wnnJ 

That  is  W"1  is  also  diagonal.  G  is  equal  to  the  multiplication  of  two 
diagonal  matrices  and  therefore  diagonal. 
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Now,  we  want  to  prove  that  if  G  is  diagonal  then  PN  is  diagonal. 

(66)  is  equivalent  to 

[  I+PNRF  ]G  =  PNR 
or  PNR  -  PNRFG  =  G 

or  PN  =  [  R-RFG  ]~lQ  (67) 

Since  G  is  diagonal,  the  right  hand  side  of  (67)  are  all  diagonal  matrices. 
Follow  the  arguments  as  the  former  prove,  we  must  have  PN  is  diagonal.  When 
the  plant  has  a  V-constra i ned  structure,  we  can  use  (56)  or  P=[I-HK]“^H  in  all 
the  proofs  and  the  same  results  hold. 

To  summarize(refer  to  Figure  26a), 
decouple  the  closed  loop  system  is 
that  the  matrix  PN  is  diagonal. 

Let  [P]q  be  the  diagonal  matrix  of  P  by  eliminating  all  the  non-diagonal 
elements,  then  a  solution  to  make  PN  diagonal  is 

N  =  P-1  [P]D  (68) 

because 

PN  =  P  P-1 [P] o  »  [P]0. 

We  are  interested  in  simple  implement  ions  of  decouplers  in  the  form  of  (68). 


with  the  matrices  R  and  F  diagonal,  to 
equivalent  to  design  a  decoupler  N  such 
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N„ 

n12  *  * 

’  N1n 

n21 

n22*  • 

*  N2n 

• 

• 

9 

• 

• 

• 

Nni 

Nnr*  * 

*  Nnn 

B. 


A. 


FIGURE  27 
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Let  us  go  back  to  Figure  27(a)  and  (b) .  We  can  reconfigure  Figure  27(a)  to 
obtain  Figure  28.  It  is  important  to  notice  that  Figure  28  is  a  feedforward 
structure  and  Figure  27(b)  is  a  feedback  structure.  For  the  simplicity  of 
implementation,  we  are  especially  interested  in  the  one  shown  in  Figure  28 
with  all  the  diagonal  elements  equal  to  unity;  and,  for  the  structure  shown  in 
Figure  27(b),  we  are  interested,  in  the  one  with  H=I.  These  decoupling 
structures  arj  shown  in  Figure  29(a)  and  (b)  respectively. 

From  Figure  29(a),  it  is  easy  to  see  that  the  P-constra i ned  decoupler  becomes 
a  feedforward  path  from  the  outputs  of  the  compensator  R  to  the  inputs  of  the 
plant  P.  From  Figure  29(b),  we  also  notice  that  the  V-constra i ned  decoupler 
becomes  a  feedback  path  from  the  inputs  of  the  plant  P  to  the  outputs  of  the 
compensator  R.  These  feedforward  and  feedback  decoupling  structures  are 
simple  implementations  of  (68).  In  the  next  section,  we  are  going  to 
characterize  these  structures  and  some  explicitly  formulae  are  also  derived. 
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4.4  The  Structure  of  Decoupled  systems 


In  the  previous  sections,  we  intended  to  decouple  a  MIMO  system  into  a  set  of 
SISO  systems.  But  what  is  the  structures  of  these  systems?  From  an 
implementation  point  of  view,  we  would  like  to  have  the  set  of  SISO  systems 
after  decoupling  equal  to  the  set  obtained  by  simply  opening  all  the  coupling 
channels  of  the  MIMO  system.  It  is  because  we  can  save  a  lots  of  time  in 
evaluating  the  resultant  set  of  systems  after  the  decoupling.  The  set  of  SISO 
systems  obtain  by  opening  the  coupling  channels  is  at  hand  well  before  the 
decoupler  is  designed.  Indeed,  we  have  virtually  separated  the  process  of 
decoupler  design  and  the  SISO  systems  design  by  this  requirement.  As 
mentioned  in  the  last  section,  we  are  only  interested  in  the  decoupler 
structures  given  in  Figure  29(a)  and  (b) .  This  section  will  give  the  explicit 
design  formulae  of  the  decouplers  and  the  limitations  are  stated. 

For  more  understanding,  consider  the  2- inputs  and  2-outputs  closed  loop  system 
given  in  Figure  30(a) 

Opening  the  coupling  channels  is  equivalent  to  make  P]_2  and  P21  zero  and  the 
equivalent  block  diagram  is  shown  in  Figure  30(b).  The  transfer  matrix  of 
this  opened  system  can  be  found  by 

Cd  =  [  I+PQRJ'^DRMd 

r  1  0  1  r  1+R11P11  0  i  r  r^^Ph  o  t 

=111  II  I  Md 

l  0  1  J  l  0  I+R22P22  -I  L  0  ^22^22-1 
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f  Rllpll 

1+R11P11 


L 


R22P22 

1+R22P22  J 


Md 


(69) 


When  a  P-constra i ned  decoupler,  as  shown  in  Figure  29(a),  is  inserted  into 
Figure  30(a),  Figure  30(c)  is  obtained.  The  transfer  matrix  is  given  by 


Cn  *  [  I+PNR  ]-lPNRMn 


S  i  nee 


r  Pn  P12  i  r  Nn  n^2  i  f  Rn  o  i 

PNR  =  |  ||  ||  I 

L  P21  P22  -I  N21  N22  J  I  0  R22J 


J  Rll(pllNll  ♦  P21N2l)  r22(P11n12  *  P12N22)  1 

L  Rn  (P21N11  ♦  P22N2l)  r22(P21n12  *  P22N22)  J 

.  [  Zl  z’  ] 
l  Z3  Z4  J 

Using  Grassmann  algebra  with 

«1T  *  t  1.0  ],  c2T  •  t  0,1  ] 

and  bj  *  [  Zj.Zj  ]T,  b2  =  [  Z2,Z4  ]T 


then 


Ad  =  det [  I ♦PNR  ] 

*  (l*Zj) (1>Z4)  -  Z2Z3 


and  we  have 
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f  (1*Z4)Z1-Z2Z3  Z2  1 

J  Ad  ld  i 

I  z3  (UZ1)Z4-Z2Z3  I 

l  Ad  Ad  J 

To  have  the  system  decoupled,  we  need 

Z2  =  Z3  =  0 

that  is 


N22p12  *  n12P11  *  0 
and  N11P2I  ♦  N2jP22  =  0 

For  our  decoupler,  Nji  =  N22  =  1,  then 


n12  =  “p12/pll»  and  N21  s  "p2l/p22- 
Now  the  decoupled  system  is 


r  R11 (PllP22"P12P2l)  1 

|  -  0  | 

p22*Rll(Pllp22-p12P2l)  I 

I  R22 (pllp22“p12p2l)  I 

l  P11*R22 (pllp22"p12p2l)  J 


(70) 


When  the  inputs  Md  =  Mn,  from  (69)  and  (70),  it  is  clear  that  Cd  t  Cn.  That 
means,  by  using  the  decoupling  structure  in  Figure  30(c),  we  cannot  make  the 
set  of  SISO  systems  described  by  (70)  equals  to  the  set  of  SISO  systems 
obtained  by  opening  all  the  coupling  channels  in  Figure  30(a). 
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The  reason  that  we  cannot  decouple  the  system  as  we  want  is  very  simple. 
Since  the  structure  of  the  decoupler  is 


r  1  -Pl2/pll  1 
1  I 
l  -P21/P22  1  J 


then 

r  pn  P12 1  r  1  -P12/P11  1 

PN  =  |  ||  | 

L  P21  p22  1  L  -p2l/p22  1  J 


r  pll“ Cp12p2l/P 22)  0  1 

I  I 

L  0  P22‘(Pl2P2l/Pll)  J 


which  is  not  equal  to  [P]q 


[P3d  = 


[Pn 

L  P2i 


p12  1 

P22  J 


Now,  examine  Figure  30(d).  We  inserted  a  V-constra i ned  structure  decoupler 
before  the  plant  matrix  P.  We  put 


r  1  0  i 

H  =  |  |  ,  and 

L  0  1  J 


r  0  -P12/P11  i 
I  l 
1  -P21/P22  0  J 


then 


PN  =  P[  I-K  ]-l 

r  pu  p12  1  r  1  P12/P1I  j*1 

L  P2i  P22  J  L  p2l/p22  1  J 
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r  pn  pi2  i  r  1  -Pi2/pn  1 
=1  II  I  (pllp22/ (pllp22“p12p2l)) 

L  P21  P22  J  l  -p2l/p22  1  -I 

r  pn  0  i 
=  I  I 

L  0  P22  J 

which  means  that  this  closed  loop  system  can  now  be  decoupled  into  the  form 
(69)  and  our  requirements  are  achieved  by  using  a  V-constra i ned  decoupler. 

Since  we  have  2  decoupler  structures  and  2  plant  structures,  there  will  be 
four  combinations  (  P-  and  P-,  P-and  V-,  V-and  P-,  V-and  V-  ) .  We  will  show 
that  in  order  to  obtain  our  desired  decoupled  systems,  we  need  to  use  a  V- 
constrained  decoupler  for  a  P-constra i ned  plant,  or  vice  verse. 

From  this  example,  we  know  that,  in  order  to  decouple  a  P-constra i ned  plant 
structure  and,  at  the  same  time,  to  fulfill  all  our  requirements,  all  we  need 
is  to  require  the  matrix  product 

pN  =  [P]D  (71a) 

but  now  N  has  only  2  allowable  form.  One  is  N j ; =1  for  all  i  and  the  other  is 
N=[I-K]~1,  where  K  has  all  diagonal  elements  zero.  To  decouple  a  V- 
constrained  plant  structure,  we  need  to  modify  (71a)  to 

[  I-HK  ] ~1HN  =  H  (71b) 

it  has  the  meaning  that  the  coupling  feedback  matrix  K  is  disconnected. 
Furthermore,  to  simplify  our  formulation,  we  need  a  short  handed  notation. 
Given  any  matrix  A,  the  matrix  [A]q  is  the  one  with  only  the  diagonal  elements 
of  A. 
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I.  P-constra i ned  decoupler  structure  vs.  P-constra i ned  plant  structure. 

The  block  diagram  is  shown  in  Figure  31.  All  we  need  is  to  find  out  whether 
(71a)  is  valid  or  not. 


We  assume  that 

P[I*N]  =  [P]D 
or  I+N  =  P-^Pb 

or  N  =  P— 1 [P] o  -  I  (72) 

but  the  right  hand  side  of  (72),  in  general,  does  not  have  all  diagonal 
elements  zero  which  is  contradictory.  For  example  take 


P  = 


r  PU  P12  1 
I  I 

L  P2i  P 22  J 


P*1[P]D  '  I 


r  Pi2P2i 

I 

1  -pllp21 


*p22p12  1 

I  (l/(pllp22-p12p2l)) 
P21P12  J 


which  does  not  necessary  have  all  diagonal  elements  zero.  So,  in  general,  we 
cannot  decouple  a  system  as  we  want  by  using  this  configuration. 
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HlaURE  31 


II.  P-constra i ned  decoupler  structure  vs.  V-constra i ned  plant  structure. 


The  block  diagram  is  shown  in  Figure  32.  From  (72a),  we  need  to  prove  that 


PC  I-N  3"1  =  [P]D 

or,  [  I-N  ]  =  [P]0-1P 

or,  N  =  I  -  [P]D-1P  (73) 


Si  nee 


[P]0 


-1  = 


r  i/pu  o  o 

j  0  1/P22  0 

I  ! 

L  0  0  ... 


0 

0 


1/Pnn  J 


and  [P]D-!P  =  | 


r  1  P12/P22 

p2l/pll  1 


1  pnl/pll 


Then  (73)  becomes 


N  = 


T  0  -Pl2/P22 

j-p2l/pll  0 

I 

l-pnl/pll 


Pln/Pnn  3 
P2n/Pnn  I 


'Pln/Pnn  3 
"P2n/Pnn  I 

I 

I 


which  means  (73)  is  always  true  and  we  can  use  this  decoupling  structure  to 
accomplish  our  requirements. 
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Pit  Pl2  • 


e 

.  .  j 


.  o 
•  o 
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FIGURE  32 


III.  V-constra i ned  decoupler  structure  vs.  P-constra i ned  plant  structure. 


The  block  diagram  is  given  in  Figure  33.  From  (72b),  we  need  to  prove  that 


[  I-HK  ]-lH[  I+N  ]  =  H 

which  is  equivalent  to 

P[  I+N  ]  =  H 
or  [  I+N  ]  =  P-!h 

and  N  =  P-1  -  I  (74) 


Since  H-*  =  [P-*]d 


It  implies  that  (74)  is  always  true  and  we  can  use  this  decoupling  structure 
to  accomplish  our  requirements. 
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•  •  0  N21  1  •  •  •  N2n  0  V22  • 


o 


FIGURE  33 


IV.  V-constra i ned  decoupler  structure  vs.  V-constra i ned  plant  structure. 


The  block  diagram  is  given  in  Figure  34.  From  (72b),  we  need  to  prove  that 

[  I-HK  ]-lH[  I-N]-l  =  H 

because  P  =  [  I-HK  ]~^H,  we  have 

[  I-N  ]-l  =  P-lH 
or  N  =  I-H-ip 


Since  H-l  =  [P"1]D  =  [  adj  P  ]D/(det  P) 


and  let  the  cofactors  be  Ajj,  then 


f  P11A11  p12a12 

,,  P12a12  p22a22 
H'lP  =  (1/detP)  j  . 

I  • 

PnlAnl  Pn2An2 


•  PlnAln  1 

•  P2nA2n  I 

I 

I 

•  PnnAnn  ^ 


As  we  known 

det  P  =  E  PijAij  =  £  PijAij 

i=l  j=l 

In  general,  we  have  (Pj jAj j)/(det  P)  t  1,  for  i=l,...,n. 

That  means,  [  I-H-1R]q  are  not  all  zeros.  It  implies  that,  in  general,  we 
cannot  decouple  a  system  as  we  want  by  using  this  configuration. 
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.  ZZA  0  0  ...  i  0  —  0  •  • 


FIGURE  34 


We  can  conclude  here  that  ,  from  our  proofs,  when  a  P-constra i ned  plant  model 
is  used,  we  need  a  V-constra i ned  decoupler  to  fulfill  all  our  requirements,  or 
vice  verse. 

As  we  mentioned,  from  only  the  inputs  and  the  outputs  of  a  MIMO  system,  there 
is  no  way  to  identify  whether  it  is  P-  or  V-  constrained.  We  need  to  assume 
that  the  structure  of  the  MIMO  system  has  been  determined  by  some  physical 
means,  and  once  this  mathematical  model  is  at  hand,  the  results  of  this 
section  can  be  applied.  Notice  that,  for  example,  when  we  have  a  P- 
constrained  plant  structure,  the  most  effective  way  to  decouple  is  to  use  a  V- 
constrained  decoupler  and  the  parameters  are  specified  by  (74).  But,  for 
certain  reasons,  we  must  use  a  P-constra i ned  decoupler  structure,  we  can 
transform  the  mathematical  model  of  the  plant  into  a  V-constra i ned  one  and  use 
this  identified  structure  to  specify  the  P-constra i ned  decoupler  parameters  by 
(73). 
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5.0  Multivariable  Dead-Time  Compensator 


In  this  section  we  will  generalize  the  well  known  Smith  compensator  to  the 
multi  variable  case. 

5.1  Smith  Compensator 


Consider  the  SISO  time-delayed  system  shown  in  Figure  35,  where  Ri(s)  and 
C^(s)  are  the  Laplace  transform  of  the  input  and  output  functions.  Tr(s), 
Tp(s),  and  Tp(s)  are  the  transfer  functions  of  the  input  compensator,  the 
plant,  and  the  feedback  compensator  respectively,  r  is  the  time  delay  of  the 
plant.  The  transfer  function  of  the  block  diagram  shown  in  Figure  35  is 


Ci(s)  TR(s)TP(s)e-rs 

Tl(s)  = - =  . - .  (75) 

Rj,  (s)  1  -  TF(s)TR(s)Tp(s)e“rs 

Let  TR(s)  =  1/ (s+3) ;  Tp(s)e-™  =  e-0-2*/[(s*l)  (s+2)] ;  TF(s)  =  1/s, 
then 

1  1 

.  e-0.2s 

(s+3)  (s+1) (s+2) 

Ti(s)  =  — . 

1  1 

1  «.  _ _ e-0.2s 

(s+3)  (s+1) (s+2) 


s  e“®-2s 

(s)  (s+1)  (s+2)  (s+3)  +  e"®  -2s 


(76) 


The  existence  of  the  time-delay  in  the  denominator 
example,  is  a  problem  of  systems  stability.  Now, 


of  Tj(s),  as  seen  from  the 
consider  the  block  diagram 
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shown  in  Figure  35  with  a  smith  compensator  added,  which  is  redrawn  in  Figure 


The  transfer  function  of  the  block  diagram  shown  in  Figure  36  is 


Tr(s) 


Tp(s)e- 


T2(s)  = 


C2(s)  1  ♦  TP(S)TR(s)TF(s)(l-e-^) 

R2(s)  '  TF(s)TR(s) 

1  ♦ . . Tp(s)e"rs 

1  *  Tp(s)TR(s)Tp(s)  (l-e~rs) 

Tr(s)  Tp(s)e_rs 


1  ♦  TF(s)TR(s)Tp(s)e- 


With  the  example’s  numerical  values, 


T2(s)  = 


(s*3)  (s+1)  (s+2) 


,-0.2s 


1  1  1 

1  . . . e-°-2s 

s  (s*3)  (s+1) (s+2) 

s  e"0-2s 

(s)  (s+1)  (s+2)  (s+3)  *  1 


Examine  the  denominator  of  T2(s),  we  find  that  there  is  no  time-delay  in  it. 
That  is,  the  problem  of  system  stability  raised  from  the  time-delay  is  solved. 
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FIGURE  36 


5.2  Multivariable  Smith  Compensator 


To  generalize  the  Smith  compensator  structure  to  the  multivariable  cases,  we 
need  to  examine  Figure  36  in  more  details.  If  we  simplify  the  Smith 
compensator  structure  shown  in  the  Figure,  we  can  have  the  simplified  block 
diagram  shown  in  Figure  37. 

The  input  and  output  of  the  Smith  compensator  is  marked  with  X(s)  and  W(s) 
respectively.  The  transfer  function 

W(s) 

SH(s)  =  -  =  (l-e_TS)  TF(s)Tp(s) 

X(s) 

=  Tp(s)Tp(s)  -  TF(s)Tp(s) e“rs  (79) 

Notice  that  the  first  term  of  (79)  is  the  product  of  the  transfer  functions  of 
the  plant  and  the  feedback  compensator  without  the  time-delay  while  the  second 
term  is  the  same  product  with  the  time-delay. 

For  the  multivariable  case,  let  the  transfer  matrices  of  the  plant  and  the 
feedback  compensator  are  P  and  F  respectively.  To  generalize  further,  we 
assume  that  there  are  different  time-delays  in  the  matrices  P  and  F.  If  the 
transfer  matrices  with  no  time-delays  (means  all  the  r* s  are  zero)  are  P*  and 
F*,  then  the  multivariable  Smith  compensator  will  have  the  structure 

SH  =  F*P*  -  FP  (80) 

notice  that  (80)  is  a  transfer  matrices  equation.  To  elaborate,  consider  the 
multivariable  system  shown  in  Figure  38.  Where  R,  P,  and  F  are  the  transfer 
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FIGURE  37 
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FIGURE  38 


matrices  of  the  input  compensator,  the  plant,  anr1  feedback  compensator 
respectively.  (s) ,  M|2(s),  •••,  Mjn(s)  ]  is  the  Laplace  transform 
of  the  input  functions  while  [  Cn(s),  C^s),  ...,  Cln(s)  ]  is  that  of 
the  output  functions.  Consider  an  numerical  example  with 


P  =  I 
I 
I 


(s+1) 


e 


■0.1s 


L  (s+3) 


(s+2) 


,-0.5s 


R  =  | 

I 

I 


1 

(s) 


F  = 


l 

r  e-°-2s 
I 

L  0 


(s+3)  J 

0  1 
e-0-3s  J 


Since  Cj  =  [  I+PRF  I'lpRMj 


r  e-o-u 


PR  =  | 


(s) (s+1) 


I  1 
I 


L  (s) (s+3) 


-0.5s 


(s+2) (s+3)  J 


f  q-O.Zs 


PRF  =  | 

I 


I  (s) (s+1) 


-0.2s 


|  - 

L  (s) (s+3) 


-0.8s 


(s+2) (s+3)  J 


(81) 
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Using  Grassmann  algebra  with 


ClT  =  [  1,0  ],  c2T  =  [  0,1  ], 


r  e-°-ls  ■ 

1 

o  I 

1  (s)(s+l) 

|  1 

,  b2  = 

e~0.5s  | 

L  (s)  (s+3) 

.  (s+2) (s+3)  J 

We  have 

r  e-°-ls  1 

I  -  o  | 

|  (s)(s+l)+e-°-3s  | 

Ci  =  |  | 

1  (s+2) [s (s+1) +e~0 • 3s_e-0 . Is]  e-0.5s  j 

L  (s) [s (s+1) +e“® ■ 3s] [ (s+2) (s+3) +e“® • 3s]  (s+2) (s+3) +e-0 • ®s  J 

(82) 


Examine  (82),  we  find  that  there  are  time-delays  in  the  denominators.  In 
general,  the  characteristic  equation  of  (81), 


det  [  I+PRF  ]  =  0 


contains  time-delays  if  matrices  P  and  F  have  time-delays.  When  multivariable 
Smith  compensator  (80)  is  used,  the  block  diagram  is  shown  in  Figure  39.  To 
illustrate,  use  the  numerical  example  again,  then  (80)  becomes 


Sp  = 


1  0 


0  1 


r  i 

1 

o  I 

I  (s+l) 

1  1 

1  | 

l(s+3) 

(s+2)J 

,-0.2s 


-0.3s 


i  r  e-°  ls  1 

.  0 

(s+1) 

1  0-0. 5s 

J  l  (s+3)  (s+2)  J 


no 


r  1  i 

|  - (1  -  e-°-3s)  0  | 

I  (s*l)  j 

|  1  1  | 

|  - (1  -  e'0-3s)  - (1  -  e-°-8s)  | 

L  (s+3)  (s+2)  J 

Simplify  the  inner  loop  of  Figure  39,  we  have  the  resultant  input  compensator 

R*  =  [  I+RSF  ]-'R 

Using  Grassmann  algebra  with 

C1T  *  t  1,0  ],  «2T  ■  t  0,1  ], 


r  i  i 

1 

1 - 1 

0  1 

1  (•)  1 

!  1  ,  b2  = 

l  1 

1  i 

1  0  1 

l  J 

(s+3)  J 

let  dj  =  (s)  (s+1)  +  (l-e“®,3s), 
and  d2  =  (s+2) (s+3)  ♦  (l-e~®-8s)  . 

Then 


R*  = 


r  (s+1)  1 

-  0 

dl 

-(s+l)(s+2)(l-e-°-3s)  (s+2) 

L  (s+3)  djd2  d2  J 
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Proceed  as  the  previous  developments  with  R*  as  the  input  compensator  instead 
of  R,  we  have 


C2  =  C  I  ♦  PR*F  3_1  PR*M2 
let  d3  =  (s+2) (s+3)  ♦  (l-e"®-^s) 

then 


r  e-O.ls 


PR*  = 


dl 

(s+1)  d3 


L  (s+3)  djd2 


1 

0  I 

e-0.5s  I 
d2  J 


T  e-^ • 3s 


PR*F  = 


dl 

(s^l)  d3 


l  (s+3)  d]d2 


1 

0  I 

e-0 . 8s  I 
d2  J 


Using  Grassmann  algebra  with 

ciT  -  [  1,0  ]#  c2T  «  [  0,1  ], 


r 

i 

e-0.1s  i 

r 

| 

1 

n  i 

i 

1 

bi  =  i 

1 

i 

dl  1 

(*♦1)  d3  | 

i 

1 

,1*2-1 

1 

i 

yj  1 

e"0-5s  | 

i 

i 

(s+3)  did?  J 

1 

l 

d2  J 

we  have 


f  e-°ls 

|  (s)(s+l)+l 

|  (s+l)[(s+2)(s+3)+l-e-°5s] 

l  (s+3)  [  s(s+l)+l  ][  (s+2)  (s+3)  +1  ] 


1 

0  I 

I 

e-OSs  | 

-  | 

(s+2) (s+3) +1  J 

(83) 


Comparing  with  (82) ,  we  find  that  there  are  no  time-delays  in  all  the 
denominators  of  (83).  This  is  because  we  use  a  multivariable  Smith 
compensator  as  shown  in  Figure  39  and  equation  (80).  A  general  prove  of  (80) 
is  developed  here  now. 


Refer  to  Figure  39  again,  but  now  ^2^=[M2i (s) , . . . ,M2n(s)]  and 
c2T=CC2l(s) » • •  »c2n(s)] •  R,  P,  and  F  are  nXn  transfer  matrices.  We  assume 
that  there  are  time-delays  in  matrices  P  and  F  only.  Let  P*  and  F*  are 
matrices  of  P  and  F  without  the  delays  respectively.  We  want  to  prove  that  by 
using  (80),  the  characteristic  equation  of  the  resultant  transfer  matrix  from 
Figure  39  consists  of  no  time-delays. 


Writing  R*  =  [  I+RSp  ]"lR  =  T“^R,  then  Figure  39  is  equivalent  to  Figure  40. 
And 


C2  =  [  I+PR*F  ]-1PR*M2 
=  [  I+PT-iRF  ]-1PT-1RM2 

Since  [  I+PT-iRF  ]  =  P[  I+T-iRFP  ]P"1 

=  PT-1[  T+RFP  JP'1 
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imp  I ies  that 


[  I+PT-iRF  ]-l  =  P[  T+RFP  J-lTP'1 

Therefore 

C2  =  [  I*PT-1RF  3-lpT-1RM2 
=  P[  T+RFP  ]~1TP-1PT-1RM2 
=  P[  T+RFP  ]-!RM2 

As  T  =  I  ♦  RSp 

=  I  ♦  R[  F*P*  -  FP  ] 

we  have  C2  =  P[  I*RF*P*  ]"lRM2 

The  characteristic  equation  of  this  multivariable  system  is 

det  [  I+RF*P*  3  =  0 

by  our  assumption,  R,  F*,  and  P*  consist  of  no  time-delays,  which  implies  that 
the  characteristic  equation  consists  no  time-delays.  And,  we  have  proved  that 
(80)  is  a  multivariable  Smith  compensator. 


117 


6.0  More  Applications  of  Grassmann  Algebra. 

In  this  chapter,  we  will  present  two  more  applications  of  the  Grassmann 
algebra.  The  algebra  is  now  used  as  an  intermediate  step  to  obtain 
simultaneous  equations  instead  of  the  symbolic  evaluations  of  determinants  as 
the  previous  chapters 

6.1  Liapunov  Equation. 

Either  in  stability  studies,  or  in  the  ouput  feedback  optimal  control,  a 
fundamental  symbolic  matrix  must  be  solved.  That  is  the  Liapunov  equation: 

AP  ♦  APT  =  -Q  (84) 

where  A  is  the  system  matrix  of  general  linear  equation 

x  =  Ax  +bu  (85) 

and  P  is  a  symmetrical  matrix  whose  parameters  are  usually  to  be  determined. 

Grassmann  algebra  is  most  helpful  in  the  investigation.  A  general  derivation 
is  shown  as  follows. 

The  frequently  used  technique  is  to  expand  (84)  into  a  set  of  simultaneous 
equations  by  using  the  Knonecker  product 

[  AT®  I  ♦  I®AT  ]  [p]  =  -[q]  (86) 
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where©  means  Knonecker  product,  the  matrices  in  detail  are  as  follows: 


r  auI  a2il 

I 

I  «12l  a22I 


[A^SI  ♦  &AT]  = 


inli  i  r  AT  0 

!  1  T 

Inol  I  I  0  AT 


.  0  1 


L  a^nI  a2nI 


ann*  -I  LOO 


If  A  is  transformed  into  the  phase  variable  form  by  using  the  Krylov 
transformation,  e.g  (87)  will  be  simplified  to 


r  0  0  ...  0  an  II  f  AT  0  .  .  .  0  1 


[aXsi  ♦  mT]  = 


1  0  ...  0  an_xI 


0  At 


10  0  ...  I  ax  I  J  100  . . .  AT  J 


and  p  =  [  PU,P12i  •  •  •  »Pln'Pl2»P22*  •  •  ■  »P2n>Pl3>  •  •  •  »Pln*P2n»  •  •  •  *Pnn  1  or  arrange 
the  columns  of  P  continuously  into  a  vector;  q  is  defined  similarly. 


Then  the  Liapunov  matrix  equation  (84)  becomes  the  following: 


W  p  =  -q 


Equation  (89)  is  nXn  algebraic  equations 
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If  Pjj  =  Pj j ,  (i^j),  or  P  is  a  symmetrical  matrix.  There  are  only 


n(n-l) 

.  (90) 

2 

independent  equations  and  the  same  number  unknowns 


W  Pr  =  -qr 


where  the  dimension  r  is  equal  to  the  shown  in  (90) 


In  deta i I ,  we  have 


r  pn  l 

I  l 

[  Qll  1 

1  Pl2  1 
j 

i  ^12 1 

1  Pl3  1  = 

1  1 

-1  <*13  1 

1  1 

1  1 

(91) 

1  i 

i  i 

Pnn  -1 

1  1 

1  1 

<■  <*nn  ■< 

How  to  find  W  by  the  computer  is  elegently  established  as  follows: 
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Assume  W  =  [  Wrs  ] ,  A  =  [  a ; j  ]  (92) 

Let  I  =  1,2, ... ,n;  J  =  I, (I+l) , . . .  ,n 

K  =  1,2, . . . ,n;  L  =  K, (K+l) , . . . ,n  (93) 

I 

I  (K,L) 

I 

(I.J)  I  wrs 


r  i=j - 

— -  W(r,s) 

=  2A(L,J) 

r  l*j 

1 

l 

i 

L  I*J . 

— -  W(r,s) 

=  A(L,J) 

I 

l 

r  k=j - 

— +  W(r,s) 

=  2A(I,I) 

l  L=J 

1 

L  K^J - 

— -  W(r,s) 

=  A(I,I)+A(J,  J) 

r  i=j - 

+  W(r,s) 

-  2A(L,J) 

r  l=j 

1 

1 

i 

l  Ij£J - 

— -  W(r,s) 

=  A(L,J) 

l 

1 

r  k=j - 

— -  W(r,s) 

=  2A(I,I) 

L  L^J 

1 

1  r  l=i  - 

— -  W(r,s) 

=  A(k,J) 

L  K*J  1 

L  L*I  - 

— -  W(r,s) 

=  0 

The  algorithm  shown  above  is  much  simpler  than  previous  works. 

When  the  Grassmann  algebra  is  applied  to  find  W  symbolically  and  A  is 
restricted  to  the  companion  form.  We  have  the  following  results  for  n=2 

r  o  ii 

A  -  I  I 

l  -32  -aj  J 
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W  =  A1®!  +  ISA1  = 


T  0  ~a2  “a2  o  1 

I  1  “al  0  "a2  I 

i  1  0  -aj  -a2  j 

10  1  1  -2ai  J 


W  in  Grassmann  terminology  is 


-2134 

+2143 

-2413 

-3142 

-4214 

♦3412 


For  n=3  where 


r  0  1  0  1 

I  I 

A  =  |  0  0  1  | 

I  ! 

I  -a3  -a2  -ai  J 


and  W  in  Grassmann  terminology  is  as  follows 
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♦  231564789 

-  231564798 

♦  231564897 

♦  231569487 

♦  231584796 

♦  231586497 

♦  231764859 

♦  231765489 

-  231765498 

♦  231769458 
+  231785496 

♦  23-786459 

-  231789456 

♦  236584197 

♦  236784159 

♦  251764893 

♦  251769483 

♦  251786493 

-  253764189 

♦  253764198 

-  253784196 

♦  256784193 

-  431569728 

♦  431569827 

-  431586729 

♦  431589726 

♦  431765829 

♦  432569187 

♦  432586197 

♦  432765189 

-  432765198 

♦  432769158 

♦  432785196 

♦  432786159 

-  432789156 

♦  436589127 

♦  436785129 

♦  451769823 

♦  452769183 

♦  452786193 

♦  453769128 

♦  453786129 

-  453789126 

♦  456789123 
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The  matrix  Riccati  equation  is  the  core  for  solving  the  optimal  problem  of  a 
linear  plant  with  respect  to  the  quadratic  performance  index.  There  are  two 
schools  of  attack:  one  is  to  consider  the  matrix  Liapunov  equation  and  based 
on  the  successive  Liapunov  equation  solutions  to  approach  the  matrix  Riccati 
equation.  The  well  known  Kleinman  technique  belongs  to  this  school.  However, 
Kleinman’s  method  is  purely  numerical  and  iterative:  every  trial  should  start 
from  the  very  beginning.  The  other  school  is  to  expand  the  matrix  Riccati 
equation  into  a  set  of  simultaneous  equations  and  keep  the  symbolic  parameters 
as  far  as  possible.  This  section  offers  an  expansion  technique  in  the  second 
school.  The  method  is  to  establish  algorithms  in  order  to  use  Grassmann 
algebra  after  the  set  of  simultaneous  equations  formed. 

consider  a  linear  plant, 

x  =  Ax  +Bu  (94) 

The  performance  index  is: 

1  r* 

J  =  —  I  (  x^QX  ♦  u^Ru  )dt  (95) 

2  JO 

Where  Q  and  R  are  symmetric  matrices. 

The  fundamental  problem  is  to  solve  the  following  Riccati  equation: 

P  =  PBR-^P  -  (PA*AtP)  -  Q  (96) 
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with  using  the  relation 


l  =  Px  (97) 

we  obtain  the  control  law 

u  =  -R-1BT^(t)  (98) 

It  is  natural  to  think  that  we  expand  (96)  into  a  set  of  simultaneous 
equations: 


p  =  Qp-otp-q  (99) 

where 

P  s  C  Pll.Pl2.Pl3* • • • >Pln*P22»P23» • • • »Pn-l,n*Pn,n  ]T  (10°) 

and 

9  =  C  911.912.913*  •••»9ln»922-P23»  - -/qn-l.n^n.n  ]T  (101) 

we  would  like  to  construct  two  algorithms  in  order  to  find  the  corresponding 
matrices  Q  and  a.  where 

a  =  [  ars  ]  r , s=l,2, . . . ,  n(n+l)/2  (102) 

Q  =  [  ars  ]  r  =1,2, . . .,n(n*l)/2 

s  =1 , 2 , . . . , n (n+1) (n^+n+2) /8  (103) 
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(I)  Formula  for  finding  matrix  a.  (APL  function:  AFA) 


a  =  f  ars  3*  r,s=l,2, . ..,n(n+l)/2 


The  values  of  ors  are  found  by  two  steps. 


(1)  form  a  matrix  V  by  writting  the  subscripts  of  the  elements  of  the  vector 
p.  (APL  function:  INDEX) 


r  i 

1  1 

f  Pn  1 

l  i 

2  1 

j 

1  Pl2  1 

1  •  1 

1  i 

n  j 

1  Pin  1 

1  2 

2  1 

P  «  1  P22  ! 

1  2 

3  1 

I  p23  | 

i  i 

|  n-1 

n-1  j 

!J*n-l,n-l  1 

|  n-1 

n  | 

i 

|Pn-l,n  1 

Let  I  =  vrl  ,  J  =  vr2 


and  K  =  vsj  ,  L  =  vs2 


(2)  Def i ne  a  function: 


r  i 

f(x,y)  =  I 

L  0 


(x=y) 

(**y) 
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Let  a  =  f(I,K)  ,  b  =  f(I,L) 

and  c  =  f(J,K)  ,  d  =  f(J,L) 

and  T  s  1  ♦  f(I,J) 

Then  the  elements  ars  is  evaluated  by  the  following  formula: 

ars  =  T{dA(K,I)+a(l-bd)A(L,  J)  +  (l-a)  (1— d)  [cA(L,I)+b(l-c)A(K,  J)]}  (104) 

To  illustrate  the  steps  shown  above,  we  give  a  simple  example  as  follows: 

If  A  is  a  2X2  matrix 

f  An  A12  1 
A  =  |  | 

L  A2j  A22  J 

then  a  is  a  3x3  matrix  because  n(n+l)/2  =  3, 
r  *11  *12  *13  1 

I  I 

a  =  |  a2i  a22  a23  | 

I  I 

L  «31  <*32  *33  J 

We  write  the  matrix  V  first: 

riii 

I  l 

V  =  I  1  2  I 

I  I 

l  2  2  J 
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I=vn=l,  J=V12=1 ,  T=l+f (I, J)=l+f (1, 1)  =1+1=2 

S=1 :  K=vii=l ,  L=vn=l 

a=b=c=d=f (1,1)=1 
l-bd=l-a=l-d=l-c=0 

substituting  the  values  of  a,  b,  c,  d  and  T  into  (104),  we  have 
an=TdA(K,I)=2A(l,l)=2An 

s=2 :  K=V2]=1 ,  L=v22=2 

a=f (I,K)=f (1,1)=1,  b=f (I, L) =f (1 , 2) =0 
c=f (J,K)=f (1,1)=1,  d=f (J,L)-f (1 ,2)=0 
l-bd=l,  l-a=0 

•  a12-Ta (1 -bd) A (L , J) =2A (2,1) =2A21 

s=3 :  K=v31=2,  L=v32=2 

a=f(I,K)=f(l,2)=0,  b=f (I , L)  =f (1 , 2)  =0 
c=f  ( J ,  K)  =f  (1 , 2)  =0 ,  d=f  ( J ,  L)  =f  (1 , 2)  =0 
a13=° 

I=v21=l,  J=v22=l,  T=l+f (I, J)=l+f (1,2) =1+0=1 

sal:  K=1 ,  L=1 

a=f (I,K)=f (1,1)=1,  b=f (I,L)=f (1,1) =1 
c=f(J,K)=f(2,l)=0,  d=f (J,L)=f (2, 1)=0 
l-a=0 

a21=Ta(l-bd)A(L,J)=(l-c)A(l,2)=A12 


s=2:  K=l,  L=2 

a=f  (I ,  K) =f (1 ( 1) =1 ,  b=f (I , L) =f (1 , 2) =0 
c=f(J,K)=f(2Jl)=0,  d=f(J,L)=f(2,2)=l 
l-a=0 

a2 2=T  [dA  (K ,  I)  ♦a  (1-bd)  A  (L ,  J)  ]  =A  (1 , 1)  +A  (2 , 2)  =An  +A22 

s=3 :  K=2,  L=2 

a=f (I , K) =f  (1 , 2)  =0,  b=f (I , L) =f (1 , 2) =0 
c=f ( J , K)  =f (2 , 2)  =1 ,  d=f (J, L)  =f (2, 2)  =1 

•  fl23=TdA(K,I)=A(2,l)=A2l 

r=3 :  I=v31=2j  J=v32=2,  T=l+f (I, J)=l+1=2 

s=l :  K=l,  L=1 

a=f  (1,^=0,  b=f  (I/L)=0 
c=f(J  ,«)=<),  d=f (J,K)=0 

/.  231=° 

s=2 :  K=1 ,  L=2 

a=f (I,K)=0,  b=f (I , L) =0 
c=f(J,K)=0,  d=f (J,L)=1 
l-a=0 

a32=TdA(K,I)=2A(l,2)=2A12 

s=3 :  K=2,  L=2 

a=f (I , K) =f (2 , 2) =1 ,  b=f (I,L)=1 
c=f (J,K)=1,  d=f (J,L)=1 
l-bd=l-a=0 

•  a33=TdA(K,I)=2A(2,2)=2A22 
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Finally  we  have 


f  2An  2A21  0  1 

I  I 

a  =  I  A12  An+A22  A21| 

I  I 

t  0  2A12  2A22J 

(II)  Formula  for  finding  matrix  fl.  (APL  function:  PRP) 

fl  =  [  wrs  ]  ,  r=l,2, . . .,n(n+l)/2 

s=l,2, . . . ,n(n+l) (n2+n+2)/8. 

The  values  of  wrs  is  found  by  three  steps: 


(1)  I  =  vrj,  J=vr2  (  the  rl^h  and  the  r2tfl  row  of  the  matrix  V) 


(2)  Construct  a  matrix  U  by  writting  the  subscripts  of  the  elements  of  vectors 
p.  (APL  function:  INDEX1) 


r 

1 

1 

1 

1  1 

r  pn 

pll 

l 

l 

i 

1 

1 

1 

2  1 

1  pn 

i 

P12 

l 

1 

1 

1 

n  | 

1  pn 

pln 

U  =  [  Ujj  ]  =  1 

1 

1 

1 

1 

2 

2  1 

p  =  1  pn 

p22 

1 

1 

2 

3  1 

1  pn 

i 

p23 

1 

1 

1 

1 

1 

n 

n  | 

1  pll 

Pnn 

1 

1 

1 

2 

1 

2  1 

1 

1  p12 

p12 

1 

1 

1 

n-1 

n 

n-1 

n-1  j 

• 

1 pn-l , n-1 

Pn-1 . n-1 

1 

n-1 

n 

n-1 

n  j 

lpn-l,n 

Pn-l,n 

l 

n 

n 

n 

n  J 

l  Pn, n 

pn,n 
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where  i  =  1,2,  .  .  ,n(n+l)  (n2«-n+2)/8; 
j  =  1,2,3,... 

Let  K=usl)  L=us2,  M=us3,  N=us4 


(3)  Define  a  function: 


r  z 

g(x,y,z)  =  I  y 
i  o 

for  example, 


(x=y) 

(x=z) 

(x^y  and  xtz) 


g(l,l,2)=2,  g (1 , 1 , 1)  =1 
g(2,l,2)=l,  g(l,2,2)=0 


Let  a=g(I,K,L),  b=g(I,M,N) 
c=g(J,K,L),  d=g(J,M,N) 
and  T=l>f (K,M)f (L,N) 


Then  «rs  can  be  determined  by 


wrs  =  [  C(a,d)>C(b,c)  ]/T 


when  lXm=Q,  C(l,m)is  defined  as  zero. 


(105) 
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An  illustrative  example: 


If  C  s  is  a  2X2  matrix 

f  Cn  C12  1 

c  =  I  I 

L  C12  C22  J 
then  Q  is  a  3X6  matrix 

r  w12  «13  w14  w15  w16  1 

I  I 

Q  =  )  «2l  w22  w23  w24  w25  w26  I 

I  i 

L  w31  w32  W33  u34  u35  11/3 g  J 

(1)  rill 

I  I 

V  =  I  1  2  I 

I  I 

l  2  2  J 

(2)  r  1  1  1  1  i 

111121 

U  =  I  112  2  1 

|12  12| 

|1  2  2  2| 

l  2  2  2  2  J 

(3) 

r=l:  I=v11=l,  J=v12=l 

s=l:  K=un=l,  L=u12=1  ,  M=u13=1,  N=u14=1 

a=g (I » K , L) =g (1 , 1 , 1) =1 ,  b=g(I,M,N)=g(l,l,l)=l 

c=g ( J » K , L) =g (1 , 1 , 1) =1 ,  d=g(J,M,N)=g(l,l,l)=l 

T=Uf(K,M)f(LJN)=Uf(lIl)f(l,l)=l>l=2 

«n=  [C(a,d)+C(b,c)]/T  =  C(l,l)  =  Cu 
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s=2 :  K=L=M=1 ,  N=u24=2 

a=9 (I , K , L) =g (1 , 1 , 1) =1 ,  b=g(I,M,N)=g(l,l,2)=2 
c=g ( J  »K , L) =g (1 , 1 , 1) =1 ,  d=g(J,M,N)=g(l,l,2)=2 
T=l+f (K,M)f(L,N)=l+f (1,1) f (1,2) =1+0=1 
«12=  [C(a,d)+C(b,c)]/T  =  C(l,2)+C (2,1)  =  2C12 

s=3 :  K=L=1 ,  M=N=2 

a=g (I , K , L) =g (1 , 1 , 1) =1 ,  b=g(I,M,N)=g(l,2,2)=0 

c=g(J,K,L)=g(l,l,l)=l,  d=g(J,M,N)=g(l,2,2)=0 

d=b=0,  C(a,d)=C(b,c)=0 

wi3=0. 

s=4:  K=u41=1,  L=u42=2,  M=u43=1,  N=u44=2 

a=g(I,K,L)=g(l,l,2)=2,  b=g(I,M,N)=g(l,l,2)=2 

c=g ( J , K , L) =g (1 , 1 , 2) =2 ,  d=g(J,M,N)=g(l,l,2)=2 
T=1 +f (K ,  M)  f  (L ,  N)  =l+f  (1 , 1)  f  (2 , 2) =1+1=2 
•  w14=  [C(a,d)+C(b,c)]/T  =  C(2,2)  =  C22 

s=S :  K=u51=1,  L=M=N=2 

a=g (I , K , L) =g (1 , 1 , 2) =2 ,  b=g(I,M,N)=g(l,2,2)=0 

d=g(J,M,N)=g(l,2,2)=0 

w15=  0 

s=6 :  K=L=M=N=2 

a=g (I , K , L) =g (1 , 2 , 2) =0 ,  b=g(I,M,N)=0 

«16=  0 
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2: 


I=v2i=l ,  J=v22=2 

sal:  K=L=M=N=1 

c=g(J,K,L)=g(2,l,l)aO,  d=g(J,M,N)=g(2,l,l)=0 
w21=  0 

s=2 :  K=L=M=1 ,  N=2 

a=g(I,K,L)=g(l,l,l)=l<  b=g(I,M,N)=g(l,l,2)=2 

c=g(J,K,L)=g(2,l,l)=0,  d=g ( J , M , N) =g (2 , 1 , 2) =1 
T=1 ♦f (K , M) f (L, N) =l+f (1 , 1) f (1 , 2) =1+0=1 
u22~  C11 

1 

s=3:  K=L=1 ,  M=N=2 

a=9 (I j  K , L) =g (1 , 1 , 1) =1 ,  b=g (I , M , N) =g (1 , 2 , 2) =0 
d=g(J,M,N)=g(2,2,2)=2 
T«l+f  (K ,  M)  f  (L ,  N)  *l+f  (1 . 2)  f  (1 , 2)  =1 
•  «23=C(a^)=C(l,2)=C12 

s=4 :  K=1 ,  L=2,  M=l,  N=2 

a=g (I >  K , L) =g (1 , 1 , 2) =2 ,  b=g ( I , M , N) =g (1 , 1 , 2) =2 
c=g ( J , K , L) =g (2 , 1 , 2) =1 ,  d=g(J,M,N)=g(2,l,2)=l 

T=l*f(K,M)f(L,N)=l+f(l,l)f(2,2)=l*l=2 
w24=  [C(a,d)+C(b,c)]/T  =  C  (2, 1)  =  C2i 

s=5 :  K=1 ,  L=M=N=2 

a=g (I » K , L) =g (1 , 1 , 2) =2 ,  b=g (I , M , N) =g (1 , 2 , 2) =0 
d=g ( J , M , N) =g (2 , 2 , 2) =2 
T=l*f(K,M)f(L,N)=Uf(l,2)f(2,2)=l 
«25=  C  (a ,  d)  =C  (2 , 2)  =C22 
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s=6:  K=L=M=N=2 

a=9  (I  i  K ,  L)  =g  (1 , 2 , 2)  =0 ,  b=g(I,M,N)=0 

W26=  0 

3:  I=v3i=2,  J=v32=2 

sal:  K=L=M=N=1 

a=9 (I , K , L) =g (2 , 1 , 1) =0, 
cag(J,K,L)=g(l,2,2)=0, 

W3i=  0 

s=2 :  K=L=M=1 ,  N=2 

3=0, 
c=0 

w32=  o 

s=3 :  K=L=1 ,  M=N=2 

a=0, 
c=0 

«33=0. 

s=4:  K=1 ,  L=2,  M=l,  N=2 

a=g (I , K , L) =g (2 , 1 , 2) =1 ,  b=g(I,M,N)=g(2,l,2)=l 

c=g(J,K,L)=g(2,l,2)=l,  d=g(J,M,N)=g(2,l,2)=l 
T=l+f (K ,M) f (L,N)=l+f (1 , 1) f (2, 2) =1+1=2 
"34=  [C (a , d)  +C (b ,  c)  ] /T  *0(1,1)  =  Cn 
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(Ill)  Function  package  for  solving  equation  (99) 


(1)  APL  function  named  RIC, 

(2)  APL  function  named  AFA, 

(3)  APL  function  named  INDEX, 

(4)  APL  function  named  IFF, 

(5)  APL  function  named  PRP, 

(6)  APL  function  named  INDEX1, 

(7)  APL  function  named  FUN, 

(8)  APL  function  named  WITH, 

(9)  APL  function  named  UPT, 

(10)  APL  function  named  PDOT, 

(11)  APL  function  named  FORMPP, 

(12)  APL  function  named  FORM. 

the  p  matrix  is  found  by  calling  the  APL  function  named  RIC. 
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7  F*-MCiPViWiQUq2',qi  iWlIiG 
r 1]  fl  SOLUTION  TO  nice ATI  EQUATIONS 
r  2 ]  ' ENTER  N - DIMENSION  OF  SYSTEM  MATRIX  A  : ' 

[3]  N*r\ 

[4]  'ENTER  MATRIX  A 

[5]  F+-AFA  A+(N,N)pti 

f  6]  ' ENTER  M-—  DIMENSION  OF  MATRIX  R  :* 

[7]  M*{.\ 

f  8 ]  'ENTER  MATRIX  R 

[9]  R*-(M,M)pU 

[10]  ' ENTER  INPUT  MATRIX  B 

f  113  W+-FRP  G<-B+.*(m)+.*bBHN,M)pO 

r  12  3  'ENTER  VECTOR  QV=QlltQY2 . Q1N.Q22  , . .  .$2  N,Qf$3 . QNN  :' 

r  13  ]  QV<n 

r 14 ]  *  ENTER  STARTING  TIME  TN  :  * 

[  15  ]  ff+Q 

T 163  ' ENTER  INITIAL  VALUES  OF  VECTOR  P  s! 

[17]  PV*-UPT(N,N)pQ 

[18]  ' ENTER  TIME  INCREMENT  ET  (  NEGATIVE  VALUES  SHOULD  BE  USED  ):» 

r  19  ]  Af+o 

r  2  O ]  ' ENTER  THE  NUM.  OF  POINTS  OF  SOLUTIONS  OF  RUNGE-KUTTA  METHOD 

[21]  NN+Q 

[22 ]  m 

[23]  L:T<-TN+ET 

[24]  Q'\<-bT*PDOT  PV+Q3*-LT*PDOT  PV+O.SxQ2+LT*PDOT  PV+O.SxQl*-LT*PDOT  PV 

[25  ]  *6)x<31+^4+2*$2+<33 

T  2  6 ]  +L  IFF  NN>I+I+ 1 
[27]  'THE  RESULT— -  MATRIX  P  IS i' 
r  2  8  ]  [J+-P«-/V  FORMP  PV 

7 


7  F+AFA  AiRiSiNiMil-.JiKiLiViBidDiEiHiT 

[I]  fl  FORMING  MATRIX  P=[a] 

12]  V*- INDEX  NHR*-l)*pA 

r 3  ]  FHM,M<- 0.5 1+ff) p  0 

[4]  10 s  r*-7TP;,9+-l] 

T 5  ]  T+h+I^J+Vt R;2  ] 

16]  Ll'.C*-J=K<-VtSitf 

r  7  ]  H*-n-C)*B*-I=L<-VF  S#] 

i8]  e<-i-b*d*-j=l 

[9]  KR-,S>T*(AlKin*D)+(AU'iJXl=K)*E)+(UlL'tn*C)+AlK'M*H)*(I*K)* 
D 

T  10]  +G1  IFF  M>S*-S+ 1 

[II]  -*r,0  IFF  M>R<-R+ 1 
7 


7  V+ INDEX  N\I\R 

[1]  V<-[2,N)p(NpI+l)  ,R<-\N 

[2]  Li  V-+-V, T 2  ](2  ,N-I)p(  (N-I)pl+I)  ,I+R 

[3]  +L  IFFiN-UzI+I+l 

r  4  ]  7+^7 

7 


7  A+B  IFF  C 
r  11  A-CpB 
7 


I 

1 


7  W+PRP  GiNiMihJiKiLiRiSiAiBiCiDiViUiVlirn 
Cl]  A  +  G  IS  A  SYMMETRICAL  MATRIX 

C 2  ]  MM*-lipV<-INDEX  NHR*-i)*pG 
[3]  NN+l*pU+INDEXl  N 
r  4  ]  WHMM.NN)  pO 

[5]  L0i'T+Vl[2  ]  WITH  lHS+l)*Vl+VlRi1 

r 6]  LUK<ri\m+viso 

C7]  A+I  FUN  K,I+Ull2l 

T  8  ]  C*-J  FUN  K,L 

C9]  /fH/l[3] 

CIO]  r-<-l+(X=W)AL=/^-£/i[4] 

r  11]  8+  r  FUN  M,N 

C  12]  -*BRl  IFF  Q*A*B*OD*-J  FUN  M,N 
C  13  ]  +/7/?0  IFF(0=A*D)*0=B*C 

[14]  -W?/?2  IFF  0=4  «<£) 

[15]  ->BR 0  WITH  WlRtttHJtAiD^iT 

[16]  BR>.t+BRQ  WITH  W^RiSl^BiCliT 

[17]  3/?l:/7[/7;5]-*-((7r4;Zl]+(jrS;(7]  )4T 

[18]  BR0:+LI  IFF  NN>S+S+ 1 

[ 19]  +L0  IFF  MM>R+R+ 1 


iMMiNNi 


T 


7 


7  U*-INDEX  1  NiItViM 
[  1]  Af*- 1 1  p  V+- INDEX  N 

C2]  {/«-( (,M,I+Q.  )pl)  ,[2]  V 

T3]  L'.U+-U,\  1](  ( J-l)  ,0)+(  (A/,2  JpTT:  J;  ] )  ,[2  ]  V 

[4]  +L  IFF  Afel+J+l 
7 

7  F«-X  FUN  YANDZ\Y\Z 
[  1]  y<-ltY4M?Z 

[-2]  PH  (Z*X=Y)+7xX=z  W*U=:jr)A X*Z+~ltYANDZ 

7 


7  4<-B  wra  (7 
r 1]  A+B 
7 

7z/prro]v 

7  VA+UPT  AihN 

[1]  VA*-\Q*NHI+-l)*pA 

[2]  L:  VA+VA,AlIi  (J-l  )  +  iiV+l-J] 

[3]  -+L  IFF  N>I+I+1 
7 


7  PD+-PD0T  PV 

[1]  PIM-iW+.xFORMPP  PV)-(F+.*PV)+QV 
7 

7  PP*-F0RHPP  PV;M1}M2  ;Q 
[  1]  Af2«-(  (  \Q)°.<\Q)xMlHQ,Q+vPV)pPV 

[2]  PP+UPT  M2*§M1 
7 

7  P«-ff  TOffllfP  P7;5;J 

[  1]  a  FORMING  SYMMETRICAL  MATRIX  P  FROM  VECTOR  PV 

[2]  FHNtN)p~l+I+l 

[3]  £:^r;(r-l)  +  i5]-HB«-ff-t-i-r)*py 

[4]  PV+B*PV 

rs3  -H  IFF  Nil+I+i 
[6|  P+-P*-(<SP)x  Qx(\N)o.-xN 

7 
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(IV)  Illustrative  example. 


For  the  gi ven  plant 

roil  r  o  i 

X  =  I  I  ♦  I  I 

l  0  -1  J  l  1  J 

and  for  the  quadratic  performance  index 

1  r*  r  2  o  i 

J  =  -  I  {  X '  I  I  X  ♦  uT(1/2) u  }dt 

2  J0  L  0  2  J 

the  input  data  are 

N=2,  Mai,  TN=0,  NN=100,  AT=-0.05 

roil  r  o  i 
A  =  I  I,  B  =  I  I  ,  R=l/2, 

l  0  -1  J  L  1  J 

r  2  0  1 

Q  =  I  I 

L  0  2  J 


The  initial  values  of  p  a  [0,0,0]. 


P-.-PTC 


ENTER  N -  DIMENSION  OF  SYSTEM  MATRIX  A  : 

I  Js 

2 

ENTER  MATRIX  A  : 

ris 

0  1  0  "1 

ENTER  M - DIMENSION  OF  MATRIX  R  s 

Us 

1 

ENTER  MATRIX  R  s 
fir 

0.5 

ENTER  INPUT  MATRIX  B  : 

(1: 

0  1 

ENTER  VECTOR  QV=QU,Q12  , . . .  ,Q1N,Q22  , . .  .Q2N ,Q33 . ,QNN  s 

f): 

2  0  2 

ENTER  STARTING  TIME  TN  : 

[]: 

0 

ENTER  INITIAL  VALUES  OF  VECTOR  P  j 
fl: 

0  0  0 

ENTER  TIME  INCREMENT  A T  (  NEGATIVE  VALUES  SHOULD  BE  USED  ) : 

fir 

"o.os 

ENTER  THE  NUM.  OF  POINTS  OF  SOLUTIONS  OF  RUNGE-KUTTA  METHOD s 
f  Is 

100 

THE  RESULT--—  MATRIX  P  IS: 

3  1 
1  l 


i 
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7,0  Conclusions 

We  have  proposed  an  effective  way  to  compute  symbolic  system  functions  using 
the  Grassmann  algebra  in  the  single  input  case  and  extend  the  technique  to  the 
multi-input  case.  The  main  idea  is  to  take  the  advantage  of  the  sparse  nature 
of  a  dynamic  system  and  the  elegant  notations  of  Grassmann  algebra. 

In  the  graph  representation,  we  followed  Coates*  sense,  that  is  the  outgoing 
and  incoming  branches  are  corresponding  the  column  vector  and  row  vector  of  a 
matrix  respectively.  However,  the  method  is  much  simpler  than  coates* 
approach. 

For  a  hybrid  system  which  contains  analog  as  well  as  digital  notes,  a  new 
formula  has  been  established.  This  technique  is  particularly  useful  in 
frequently  encountered  computer-controlled  systems. 

A  fundamental  technique  in  multivariable  systems  analysis  and  design  is 
decoupling.  The  P-constra i ned  and  V-constra i ned  decoupling  structures  are 
studied  in  detail  under  the  light  shed  by  the  Grassmann  algebra. 

The  typical  and  well  known  Smith  pre-estimator  and  controller  is  investigated, 
particularly  in  the  multivariable  case  by  use  of  idea  of  Grassmann. 

Finally,  for  showing  the  power  of  Grassmann,  we  applied  to  the  solutions  of 
two  basic  equations.  One  is  the  Liapunov  matrix  equation  for  stability 
studies  and  one  is  the  Riccati  matrix  equations  for  optimal  control. 


142 


These  appendices  present  two  other  methods  for  symbolic  evaluations.  By 
comparing  them  with  the  Crassmann  Algebra  method,  we  can  conclude  that  the 
Grassmann  method  is  more  elegant  and  easier  to  use. 
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To  accomplish  symbolic  evaluations,  we  can  use  methods  other  than  Grassmann 
algebra.  This  chapter  is  to  present  a  symbolic  evaluation  method  based  on  the 
Fourier  transforms.  By  comparing  the  Fourier  transform  methodology  with  that 
of  the  Grassmann  approach,  it  finds  that  the  later  one  is  more  simple  in 
manipulations  and  applicable  to  evaluation  other  than  transfer  functions. 

Al.l  Symbolic  Evaluation  Via  Fourier  Transform. 

Lee[16]  developed  a  Fourier  method  to  determine  the  characteristic  function 
f (s)  =  det[  sI-A  ]  a  sN  ♦  ctN-isN"1  ♦  •••  ♦  ajs  ♦  oq  (106) 
from  a  state  space  model 

x  =  Ax  (107) 

The  method  can  be  summarized  as  follows: 

Let  W  =  exp  (2x/N) j  be  the  Nth  root  of  unity,  then  from  (106)  we  obtain 
f(Wm)  =  det[  W"mI-A  ]  (108) 

=  W"mN  ♦  aN_iWm(N_1)  ♦  ...  +  aiW""*  ♦  a©  (109) 

m  =  0,1, .. . ,N-1 


14A 


Since  m  and  N  are  integers,  therefore 


W-mN  .  (wN)-m  =  exp  (_2mT)j  =  1 

From  (109) ,  we  have 

det[  WmI  -  A  ]  -  1  =  ♦  ...  ♦  a1W~m  +  a0  (110) 

Let 

am  =  det [  W-»I  -  A  ]  -  1  (m) 

from  (110)  and  (111),  we  have 

N-l 

am  =  E  a|<W_mk  for  ms  0,1,2, N-l  (112) 

k=0 

Therefore 

N 

ak  < . >  am 

forms  a  discrete  Fourier  series  pair 
The  inversion  formula  gives  that 
1  N-l 

ak  =  —  E  amWmt<  for  ks0,l,2, .  .  .,N-1  (113) 

N  m=0 

Now  we  can  extend  Lee  s  Fourier  method  to  transfer  function  evaluation  by  use 
of  (14) 
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1 1  I ustrat i ve  examp  I e : 


For  the  given  system 

r  -1  -2  -2  1  f  2  1 

.1  III 

x  =  j  0  -1  1  |  x  ♦  j  0  j  u 

I  111 

LI  0  -1  J  l  1  J 

y  =  [  i,  i,  o  ]  x 

its  corresponding  transfer  function  is  desired. 

We  first  evaluate  the  characteristic  equation 
f(s)  =  det [  sI-A  ] 

r  s*l  2  2  1 

I  I 

=  |  0  s+ 1  -1  | 

I  I 

L  -1  0  s«-l  J 

Let  N=3,  then  W  =  e  J (2t/3) 
substituting  into  (111)  gives 
am  =  det[  W-mI-A  ]  -  1 

f  e-j (2mr/3)  ♦  i  2  21 

=  |  0  e~J (2mx/3)  +  i  _1  |  -  l 

I  -1  0  e“J (2mr/3)  ^  j  J 


146 


Then  we  have 


ao  =  13 

aa  =  1  -  j 1.732 
and  a2  =  1  ♦  j 1 .732 


using  the  inverse  formula 


N-l 

E  amW»l< 


m=0 


(113) ,  we  have 
1  2 

=  -  l  am  eJ (2mkT/3) 

3  m=0 


(k=0,l,2) 


therefore 

ag  =  5 
a^  =  5 
a2  =  3 


The  characteristic  equation  is  obtained  as  follows 

f(s)  =  ♦  3s2  ^  5s  +  5  (114) 

Then  we  try  to  find  the  function: 
f (s)  =  det[  (sI-A)  +  bc^  3 

f  s+3  4  21 

I  I 

=  |  0  s+1  -1  | 

I  I 

10  1  s+1  J 

It  is  still  a  determinant  evaluation  problem. 

Again  let  N=3  and  W  =  e  j (2*/3) 
substituting  into  (111)  gives 
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am  =  det[  W“mI-A  ]  -  1 


r  e“ J (2mx/3)  +  3  4  21 

=  |  0  e_i (2mf/3)  +  1  -1  |  -  1 

[  0  1  e-j(2mr/3)  ♦  1  J 

After  evaluation: 

3q  =  19 

ai  =  -0.5  -  j2.598 
and  32  =  -0.5  *  j2.598 

using  the  inverse  formula  (113),  we  have 

1  N-l  1  2 

ak  = -  Z  amWmk  = -  I  am  eJ (2mkr/3)  (k=0,l,2) 

N  m=0  3  m=0 

therefore 

*0  =  6 
fflj  =  8 

a  2  =  5 

Therefore  det[  (sI-A)  +bc^  J  =  s3  ♦  5s2  ♦  8s  ♦  6  (115) 

substituting  (114)  and  (118)  into  (14)  gives 

s3  +  5s2  +  8s  ♦  6 

c 1  (si -A)  - 1  b  . . .  -1 

s3  *  3s2  ♦  5s  +  5 

2s2  ♦  3s  ♦  1 
s3  ♦  3s2  ♦  5s  ♦  5 

which  is  our  desired  transfer  function. 
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A  program  based  on  the  discrete  Fourier  method  to  evaluate  the  corresponding 
transfer  function  of  a  state  space  model  is  written.  The  language  used  is 
APL;  the  program  involves  symbolic  determinant  evaluation  only.  All  the  non¬ 
zero  elements  in  the  two  determinants  are  complex  in  general.  The  program  and 
outputs  for  the  illustrative  example  are  given  in  the  next  pages 
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VTRANSFERl  □]  V 
7  TRANSFER 

[1]  n  EVALUATING  TRANSFER  FUNCTIONS  FROM  STATE  SPACE  MODELS 
121  n  VIA  FOURIER  TRANSFORM 

[3]  n  FROM  DX=AX*BU  AND  Y=CTX  TO  Y(S)/U(S)  DETERMINATION 

[4]  ' STEP  1.  EVALUATE  THE  COEFS.  OF  DET(SI-A)' 

1 5]  '  ' 

r  6]  VC+CHARA  A 
r7]  '  • 

r 8]  'STEP  2.  EVALUATE  THE  COEFS.  OF  DET(SI-A+BCT) ' 
f 93  '  ' 

C  10 J  NI+-CHARA  A-B+.*$C 

nu  •  ' 

f  121  'THE  COEFS.  OF  NUMERATOR  ARE  :  '\wNl-DC 

r 13 J  'THE  COEFS.  OF  DENOMINATOR  ARE:  ' ,fDC 

7 

NCHARAlUW  ■ 

7  CF*-CHARA  A\N\I\ X\  XI ;  50 ;  SI\ M\  V\BM\AK ;  U\J\ CA 

[I]  .Y«~(o2)itfntp/l  . 

F  2]  BM*-AK+(M+2*N)pI<-Q 

[3]  CA+(  V+Mp  1  0  )\>1 

[4]  LiUH  2*1)1 \ 2*1+1 

[5]  ,S0<-(  2oXl ) , loYl-*-JxY 

[6]  Sl*-(N,M)pSO,MpO 

[7]  BMLUlH’l  0)+CDET  SI-CA 

T 8]  -*■(  (N-l)>I+I+l)/L 

T 9]  */I(0 )  ,/I(  1 ) . A{N-1)  :  (  COMPLEX  NUMBERS  )' 

[10]  'REAL  PART  IMAG.  PART' 

111]  t>(JV,2)pBM 
[  1?.]  I<0 

r  13]  7,1  :,/<-() 
r 14 J  V*  (  2*1) + 1 2*1*1 

f  15]  L2iAKlU]<-AKlU]+mi(2*J)l\2*J+\]  CMP  SOH  2oXl )  ,loXl*--X*l*J 

r i6 ]  -»•< {N~\)>j<rjn)/L2 

r  17  ]  -*(  {N-l)>I<rI+l)/Li  . 

[  18 ]  '  ' 

fl9]  '  THE  COEFS.  ARE  IN  DESCENDING  ORDER  : '  ,fCF*-<t?(  V/AKiN)  ,1 
7 

V(7/;£T[D]7 

7  Z+-CDET  A\B\P\I\U\N',K',AxV 
[  1  ]  tM"  1  0)+2*M 

r  2]  7.< -  1  0 

[3]  L:  AT*-' MODULUS  '  CREPE  AT,  Al  ;U] 
r  4]  r=r*-Ai\r/Ai)/LL 
f  5]  Alf.P;  ]+A[P,I;] 

IS]  7,*-Z 

I” 7 3  M,iZ+Z  CMP  B*-AlhU ] 
r  8  J  ->(  (l  =  x/Z=0)vl=^lt/f-pi4)/0 
T9]  ->(1=*/B=  0  0  ) /BR1 

[10]  B<-(  1  -1  )*£?i+/B*2 

[II]  BRUAA+'B  CMP'  CREPE  AT,  Al \U] 

F 12]  A*-Vp~l*K+l 

ri.l]  LR : AlK 2* (  2*K-1  )*AA)  CMP  '  C REPEAT  All  \] 

114]  +USkK+K+\)/lB 

[15]  A*-  1  2  IA-A 

[16]  +L 

T 17]  nB7/1£{//13’B5  /I  COMPLEX  DETERMINANT 
7 
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MlODULUSlVlV 
7  RH10DULUS  7 
I'll  /M  +/Z*2)*0 . 5 

n 

i 

1  CREPE AT[\J\V 

7  CV+FUNCNAME  CREPEAT  ARGVEC 

(  I  I  CT*-,o 

12)  ->(o=P/i/?(;^r)/o 

r  3  J  CV<-tFUNCNAME,  '  2M/?Cra7’ 

fU]  C\'+CV .FUNCNAME  CREPEAT  2* ARGVEC 

V 

76’WP[m7 

V  Z^P  CT/P  Q 

II]  Z<-(-/P*<?)  1+/P’"tv? 
r  2  ]  nP  M0  $  ARE  COMPLEX  NUMBERS 
7 


l  I  1  u:;t. rati ve  Example 


,1<-3  3p~l  '2  "2  0  "1  1  1  0  "1 
0<-3  lp2  0  1 
f>3  Ipl  1  0 
TRANSFER 

RTF’’  1.  EVALUATE  THE  COEFS.  OF  DET(SI-A) 

A  ( 0 )  ( 1 ) . /I  ( 2V-1 )  :  (  COMPLEX  NUMBERS  ) 

REAL  PART  IMAG .  PART 
13  0 

1  ~1.73  2 

1  1.73  2 

THE  COEFS.  ARE  IN  DESCENDING  ORDER  :1  3  5  5 

STEP  2.  EVALUATE  THE  COEFS.  OF  DETiSI- A+BCT) 

/1(D). /1(1) . A(N-l)  :  (  COMPLEX  NUMBERS  ) 

REM,  PART  IMAG.  PART 

n  0 

"0.5  ’2.598 

’0.5  2.598 

THE  COEFS.  ARE  IN  DESCENDING  ORDER  :1  5  8  6 

THE  COEFS.  OF  NUMERATOR  ARE  :  0231 

THE  COEFS.  OF  DENOMINATOR  ARE:  1355 
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A2.0  Symbolic  Function  Generation  Using  Number  Theoretic  Transform. 


A2.1  Introduction. 


This  appendix  shows  that  orthogonal  functions  which  can  be  generated  by  a 
cyclic  group  under  multiplication  are  suitable  for  the  generation  of  symbolic 
network  functions  in  computer-aided  analysis  and  design  of  electronic 
circuits.  A  multi -dimensional  Number  Theoretic  Transform  is  defined  which 
gives  an  efficient  method  to  generate  symbolic  network  functions  without 
quantization  error. 

a  network  function  is  usually  obtained  from  a  physical  model  of  the  electronic 
circuit  concerned.  In  general,  the  function  may  be  represented  by  a  ratio  of 
polynomials  in  L,  C,  the  complex  frequency  s  and  other  parameters.  These 
polynomials  are  the  expansions  of  two  determinants  which  are  found  by  nodal, 
loop  or  other  kinds  of  network  analysis.  This  shows  that  the  problem  of 
generating  symbolic  network  functions  of  a  given  network  reduces  to  the 
problem  of  generating  the  expanded  result  of  given  determinants. 


Y(x1,x2, . . . ,xq) 


Ni-1  Np-l  Na-1 
=  E  E  ...  E 
nj=0  n2=0  nq=0 


cnl , n2, . . . , nqxlnlx2n2 • • • xqnq 


(116) 


where  xi,x2,...,xq  are  symbolic  entries  of  L,C,s,... 

and  N^-l ,N2-1 , . . .Nq— 1  are  maximum  possible  nj,n2,..,nq  respectively. 
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A2.2  Orthogonal  Functions. 


We  assume  that  N  functions  ^k.o(mk)»  ^k,l(mk)>  ^k,N-l(mk)  are  mutually 

orthogonal  between  0  and  (N— 1) .  The  function 


Nrl  N2-l  Ng-1 

Y  (1*1 ,  i»2  *  •  •  • » IBq)  =  ^  E  ...  L  cni  j  n2> .  .  .  t  nq^l,  nl  (ml)  ^2,  n2(m2)  •  •  •  ^q,  nq(mq) 

n^=0  n2=0  nq=0 

(117) 

can  be  represented  by  choosing  the  c’s  according  to  the  following  equation: 


Ni-1  No-1  Ng-1 

cnl , n2, . . . , nq  ~  *  l  ...  l  Y(mi,m2, . . ,mq)^i  1  02,n2  *  •••  ^q,nq  * 

m^=0  m2=0  rriq=0  [ 

(118) 

Equation  (116)  is  a  general  representation  of  the  symbolic  network  function  of 
a  given  network.  In  order  to  make  equation  (117)  to  be  compatible  with 
equation  (116),  let  xknk  be  represented  by  ^k,nk-  Since  the  set 
{xknk:nk=0,l, . . .N-l>  can  be  generated  by  xk°,  xk°.xk,  xk.xk,  xk2.xk, 
xk3.xk,  .  . . .  The  orthogonal  function  ^k.nk  sets  must  also  be  generated  by  a 
cyclic  group  under  multiplication,  say  ak’s  f°r  example.  That  is 


(“k)n  =  ^k,nk 

This  forms  the  basic  criterion  used  in  this  paper  for  employing  orthogonal 
functions  to  find  symbolic  network  functions  to  find  symbolic  network 
functions.  The  coefficients  «nl , n2,  , nq  'n  equation  (116)  is  equivalent  to 

the  corresponding  coefficients  in  equation  (117)  if  xknk  is  replaced  by 
lk,nk(mk)-  All  combinations  of  Y  (m^  ,m2, . . .  ,n»q)  could  be  calculated  by 
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substituting  xn^n  by  ^k,nk(mk)  *>th  all  possible  combinations  of  mj ,n»2, . . . ,mq 
in  the  determinant  under  question.  Hence  the  coefficients  ( n2, . . . nq >s 
could  be  found  from  equation  (118). 

For  the  complex  number  field,  an  =  e“j2T/N  js  the  generator  for  a  set  of 
orthogonal  function  used  to  generate  symbolic  network  functions.  This 
formulation  is  equivalent  to  doing  the  calculation  by  multi-dimensional 
Discrete  Fourier  transforms (DFT) .  In  a  finite  ring,  we  may  also  find  a  cyclic 
group  under  multiplication  modulo  a  prime  number  to  form  an  Number  Theoretic 
Transform  which  may  be  used  to  generate  the  network  function.  If  we  use  the 
DFT  to  generate  the  symbolic  network  functions,  calculation  error  can  hardly 
be  avoided,  since  both  the  DFT  and  the  Inverse  DFT  involve  complicated 
multiplications  of  irrational,  complex  numbers.  In  next  section,  we  discuss 
to  use  the  Number  Theoretic  Transform  to  generate  symbolic  functions.  This 
new  method  gives  an  error-free  calculation  of  symbolic  functions. 

A2.3  Generation  Via  Number  Theoretic  Transf orm(NI I ) . 

Let  M  be  the  modulo  base  of  the  NTT  to  be  used.  A  multi-dimensional  NTT  pair 
can  now  be  defined  as 


Nrl  N2-l  Nq-1 

Y(mj,m2, . . .  ,mq)=<  E  E  ...E  a(n^,n2,  • .  .nq)am^n^  ...  am<ln<l  >|y| 

nj=0  n2=0  nq=0 

(119) 

and 
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(120) 


Na-1  No-l  Ng-1 

■(nl>n2» • • nq)=<  Ni_1N2”^  ...  Nq-1  I  £  . .  .  2 

nj=0  n2=0  Oq=0 


a  -mlnla-m2n2  .  .a-mqnq>^ 


where  ni,n2, . . . ,nq=0,l, . . . # (Nj-1) , . . . , (Nq-1)  respectively. 


To  demonstration  the  application  of  the  NTT  to  symbolic  network  function 
generation,  let  us  find  the  input  impedance  of  the  circuit  shown  in  Figure  41. 
The  mesh  matrix  of  the  circuit  can  be  written  as 


Figure  41 


r  V;  1 

1  1 

f  l(Kr 

-10 

-r  1 

I 

r  Ij  1 

I  I 

!  1 

1  0  1  = 

1  1 

| 

1  -10(l*gm) 
i 

13+10gffl 

1 

-1  1 
j 

l  12  1 

1  1 

1  1 
l  0  J 

l  10gm-r 

-(l*10gm) 

l*21r  J 

1  1 
t  I3  J 

(121) 

where  r  and  gm  are  network  symbols. 


Let  Ms65537,  Nj=N2=4  and  the  root  of  unity,  a,  be  2®.  Equation  (119)  and 
(120)  become 
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3  3 

Y("»l,n>2)=<  E  E  a (nl > n2)28 (mlnl+m2n2) >65537 

n^=0  n2=0 


3  3 

a(ni,n2)=<4“l4“l  E  E  Y  (m^  ,m2)2“8  (m^n^+n*^n^)  >65537 
m^=0  m2=0 


substitute  r  by  28n,l  and  gm  by  28m2  into  equation  (121), 


(122) 


(123) 


A(r,gra)  =>  A(28">1,28<"2)  = 


r  io*28ml 


_28ml  1 


<  j  -10(l+28m2)  13+10. 28m2  -1  j 

L  10 . 28m2-28ml  -(l+10.28m2)  1+21. 28ml  J 


>65537 


(124) 


for  mi,m2=0,l,2,3  respectively. 

Equation  (124)  is  actually  a  physical  representation  of  the  general  equation 
(119)  or  equation  (122),  hence 

A(28ml,28m2)  =  Y(mi,m2)  and 

{Y(m1,m2) :m1,m2=0,l,2,3>  =  {  1122,57222,682,10119,32838,42235,22998,13601, 

65375 , 45738 , 6501 5,19115,31819,51496, 42459 , 22782  > 


The  inverse  transform  can  be  found  by  equation  (123)  and  this  gives  the 
coefficients  of  the  expansion  of  equation  (121). 


{  a(n!,n2) :n1,n2=0, 1,2,3}  =  {  20,0,0,0,622,20,0,0,260,200,0,0,0,0,0,0} 
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Hence  A^g,^  =  20+622r+260r2*20rgm+200r2g|n 

Similarly,  An  can  be  found  by  the  same  procedure.  We  obtain 


V;  A  20*622r*260r2*20rgm*200r2gm 
10  Ii  An  12+273r+210rgm 

A2.4  Conclusion. 

The  number  Theoretic  Transform  gives  a  completely  error  free  calculation  of 
symbolic  network  functions.  The  major  difficulty  of  the  NTT  is  that  the 
dynamic  range  is  limited  by  the  base  of  the  module  arithmetic  involved.  This 
problem  can  be  overcome  by  using  a  large  modulus  as  the  base  of  using  the 
Chinese  Remainder  Theorem. 
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